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ABSTRACT 


The  potential  of  a  relatively  new  method  to  remotely  measure  near-surface  cur¬ 
rents  and  dominant  wave  spectra  using  Interferometric  Synthetic  Aperture  Radar 
(INSAR)  is  demonstrated.  INSAR  consists  of  a  single  conventional  SAR  augmented 
by  an  additional  receiving  antenna.  The  phase  difference  between  the  correspond¬ 
ing  SAR  image  scenes  observed  by  the  antennas  provides  an  interferogram  that  is 
directly  proportional  to  the  ocean  surface  velocity  field.  This  direct  motion  de¬ 
tection  by  INSAR  suggests  a  significant  advance  compared  with  conventional  SAR 
where  the  response  to  the  moving  ocean  surface  is  indirectly  related  to  the  complex 
modulation  of  the  surface  reflectivity  by  longer  waves  and  currents. 

An  experiment  using  an  airborne  implementation  of  INSAR  to  measure  ocean 
surface  currents  and  wave  fields,  compared  with  simultaneous  ground  truth  measure¬ 
ments  using  Lagrangian  drifters  and  wave  array  data  was  conducted  in  Monterey 
Bay.  INSAR  measured  mean  current  magnitude  estimates  agree  to  within  10  per¬ 
cent  compared  with  conventional  measurements.  The  INSAR  image  wavenumber 
spectrum  is  consistent  with  the  in  situ  directional  spectrum  and  with  predicted 
numerical  refraction  model  outputs.  The  wavelength  of  the  observed  swells  are  in 
better  agreement  (correlation  better  than  0.9)  than  wave  direction. 

An  attempt  to  estimate  the  scene  coherence  time  for  L-band  SAR  was  made  by 
taking  advantage  of  the  almost  simultaneously  acquired  SAR  and  INSAR  images. 
The  obtained  mean  scene  coherence  time  0  (100  msec)  is  consistent  with  sparse 
observed  estimates  in  the  literature.  This  limited  radar  temporal  coherence,  caused 
by  the  velocity  spread  of  short  waves,  degrades  the  azimuthal  resolution  of  SAR 
and  INSAR  and  depends  on  the  sea  state  and  radar  wavelength.  The  experimental 
results  show  that  the  finite  scene  coherence  time  has  a  dominant  role  on  the  distor¬ 
tion  of  the  INSAR  (and  SAR)  image  spectrum  relative  to  the  ocean  wave  spectrum. 
The  present  study  introduces  limitations  of  interferometric  SAR  configuration  in 
imaging  nonstationary  scenes  like  the  ocean  surface. 
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I.  INTRODUCTION 


A.  INTERFEROMETRIC  RADAR  REMOTE  SENSING 
OF  OCEAN  SURFACE  VELOCITY  FIELDS 

Historically,  oceanography  has  depended  on  sparse  observations  made  from 
ships.  Now,  orbiting  platforms  offer  global  coverage  of  the  world’s  oceans  using 
various  types  of  sensors  in  a  broad  spectrum  (visible,  IR,  microwave).  “How  relevant 
are  these  surface  observations  in  understanding  today’s  oceanographic  problems? 
How  well  matched  to  the  oceanic  time  and  space  scales  are  the  satellite  samples? 
What  oceanic  processes  can  be  observed,  and  in  what  way?  And,  finally,  what 
instruments  and  systems  exist  to  provide  oceanic  data  from  space?”  (Stewart,  1985). 

Some  of  the  answers  to  these  questions  are  addressed  in  this  study  by  ap¬ 
plying  a  new  method  of  remotely  measuring  ocean  surface  velocity  fields,  using  a 
conventional  Synthetic  Aperture  Radar  (SAR)  augmented  with  two  spatially  sep¬ 
arated  antennas  (Goldstein  and  Zebker,  1987).  By  incorporating  the  physics  of 
interferometry,  this  unique  remote  sensing  method  offers  an  opportunity  for  esti¬ 
mating  synoptic  spatial  structures  of  ocean  near  surface  currents  and  wavefields  with 
high  grid  nominal  resolution  (12  m  in  this  work)  and  large  coverage,  never  before 
obtained. 

The  objectives  of  this  research  are  to  develop  a  theoretical  foundation  for 
Interferometric  SAR  (INSAR)  in  imaging  the  free  surface,  distributed  moving  ocean 
surface,  to  verify  its  feasibility  with  simultaneous  ground  truth  experiments,  and  to 
develop  interferometric  SAR  image  processing  techniques  to  estimate  spatial  synop¬ 
tic  currents  and  wave  spectra.  An  experiment  was  conducted  to  ground  truth  the 
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velocity  field  imagery  of  INSAR  using  surface  Lagrangian  drifters,  wind  measure¬ 
ments,  and  a  shallow  water  wave  array. 

B.  STATEMENT  OF  OBJECTIVES  AND  PROBLEMS 

Three  specific  objectives  are  outlined  for  this  dissertation: 

1.  To  develop  a  theoretical  basis  for  the  capability  of  SAR  in  an  interfero¬ 
metric  mode  to  measure  the  ocean  surface  velocity  field.  The  basic  limitations  in 
correlating  the  backscatter  radar  echo  with  ocean  phenomena  is  that  the  ocean  sur¬ 
face  combines  random  fluctuating,  non-station  ary  and  non-homogeneous  scattering 
elements.  In  this  study,  a  preliminary  theoretical  model  is  developed  which  repre¬ 
sents  the  major  physical  parameters  involved  in  remotely  measuring  near  surface 
currents  and  wavefields  using  synthetic  aperture  radar  in  an  interferometric  mode. 

The  decorrelation  of  the  ocean  surface  scene  (finite  scene  coherence  time) 
as  a  function  of  INSAR  parameters  is  quantified.  This  allows  to  characterize  the 
sensitivity  and  capability  of  INSAR  in  terms  of  imaging  near  surface  ocean  currents 
and  wave  fields.  Noncoherent  motion  effects  on  airborne  SAR  (and  hence,  INSAR) 
imagery,  which  is  an  area  of  controversy,  is  also  addressed  in  the  present  study. 

2.  To  develop  analysis  techniques  to  obtain  quantitative  spatial  information 
of  near  ocean  surface  currents  and  wave  fields  from  INSAR  imagery.  Unlike  SAR, 
which  has  serious  shortcomings  in  obtaining  directional  wave  spectra,  INSAR  has 
an  inherent  improved  potential  of  measuring  2D  wavenumber  spectrum  of  the  ocean 
surface  dominant  wave  field.  This  improvement  is  related  to  the  direct  motion  de¬ 
tection  mechanism  of  INSAR  when  sufficient  Bragg  scatterers  are  present.  In  this 
study,  appropriate  theoretical  tools  are  developed  to  separate  wave  field  velocities 
from  ocean  currents.  This  separation  allows  accurate  and  unique  information  con¬ 
cerning  spatial  structures  of  dominant  ocean  wavefields,  currents,  and  other  features. 
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3.  To  ground  truth  the  currents  and  wave  fields  measured  by  INSAR.  A  si¬ 
multaneous  in  situ  and  radar  remote  sensing  experiment  was  conducted  to  determine 
the  feasibility  of  INSAR  to  image  near  surface  ocean  dynamics  compared  with  more 
conventional  measurements.  The  analyzed  data  also  provide  a  unique  study  of  wave 
refraction  and  shoaling  and  currents  in  the  complex  nearshore  region  of  Marina  in 
Monterey  Bay. 

C.  SUMMARY 

The  present  study  demonstrates  the  feasibility  of  this  unique  technique  to 
provide  reliable  quantitative  spatial  information  about  dynamic  processes  at  the 
ocean  surface.  Realization  of  the  potential  of  INSAR  with  spacecraft  might  make 
it  possible  to  routinely  monitor  global  ocean  wave  spectra  and  surface  currents.  A 
review  of  electromagnetic  interaction  with  oceanic  backscattering  waves  is  presented 
in  Chapter  II  along  with  backscattering  models  which  supports  the  interpretation 
of  the  results.  The  principal  mechanisms  of  SAR  and  INSAR  and  relevant  theory 
in  imaging  ocean  surface  dynamics  are  presented  in  Chapter  III.  In  Chapter  IV,  a 
preliminary  theoretical  analysis  for  INSAR  in  imaging  ocean  surface  velocity  fields  is 
presented  for  the  first  time.  Based  on  the  SAR  two-scale  model,  a  theoretical  model 
for  estimating  the  temporal  decorrelation  of  the  ocean  surface  using  simultaneous 
SAR  and  INSAR  observations  is  also  developed.  Chapter  V  describes  oceanographic 
relevant  theory  which  is  essential  for  the  interpretation  of  the  results.  A  numerical 
refraction  model  (Dobson,  1967)  is  used  to  predict  wave  length  and  direction  in 
shallow  water.  These  outputs  are  compared  with  INSAR  spectral  analysis  estimates 
(Chapter  VII).  Two-dimensional  wavenumber  spectra  for  in  situ  and  interferometric 
SAR  are  successively  described.  A  description  of  the  experiment,  environmental 
conditions,  and  data  processing  follows  in  Chapter  VI.  Data  analysis  methods  and 
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experimental  results  are  presented  in  Chapter  VII.  Results,  analysis,  and  discussions 
are  given  in  Chapter  VIII.  These  findings  are  summerized  and  the  major  conclusions 
and  suggestions  for  further  study  are  stated  in  Chapter  IX. 

J 
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II.  LITERATURE  REVIEW 


A.  INTRODUCTION 

Due  to  rapid  advances  in  microwave  remote  sensing  of  ocean  surface  during 
the  last  decade,  much  progress  has  been  made  in  understanding  the  physics  of 
scattering  of  microwave  radiation  from  the  ocean  surface  (Hasselmann  et  al.,  1985; 
Thompson,  1989).  The  result  of  this  progress  yields  a  capability  to  observe  in  more 
detail  near-surface  ocean  features  such  as:  wave  height,  wind  speed,  spatial  wave 
field  data,  tides,  eddies  and  currents  (Valenzuela,  1989). 

The  interpretation  of  microwave  remote  sensing  data  is  done  through  the 
physics  of  microwave  scat  tering  theories.  Therefore,  the  discussion  in  this  chapter 
begins  with  the  mechanisms  of  electromagnetic  scattering  from  ocean  surface  for 
real  aperture  radars,  which  is  the  basis  of  radar  scattering  theory.  Next,  the  dis¬ 
cussion  extends  to  coherent  synthetic  aperture  imaging  radars,  which  in  turn  adds 
the  comple  motion  effects  when  sensing  the  ocean.  The  complexity  is  mainly  a 
result  of  the  random  uncorrelated  motion  of  individual  scatterers,  smaller  than  the 
spatial  resolution  cell  of  SAR  and  of  the  order  of  the  radar  wavelength.  The  inter¬ 
action  between  different  scales  of  waves  and  currents  smear  and  defocus  the  SAR 
image,  resulting  in  difficulties  in  correlating  the  image  with  the  sensed  phenomena 
(see  Chapter  IV).  In  spite  of  the  difficulties,  considerable  progress  has  been  made  in 
understanding  SAR  imagery  of  ocean  surface  features.  (Raney,  1971;  Elachi,  1977; 
Alpers  and  Rufenach,  1979;  Valenzuela,  1980;  Plant  and  Keller,  1983;  Tucker,  1984; 
Hasselmann  et  al.,  1985;  Kasilingham  and  Shemdin,  1988). 
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B.  ELECTROMAGNETIC  FIELDS  SCATTERED 
FROM  THE  OCEAN  SURFACE 

The  interaction  of  electromagnetic  waves  with  ocean  surface  is  a  complex 
boundary  value  problem,  for  which  exact  solutions  are  not  available  to  the  present 
time  (Barrick,  1968;  Valenzuela,  1971).  This  is  because  the  ocean  surface  is  a  com¬ 
plex  distributed  scatterer,  varying  in  time  and  space,  and  also  weakly  nonlinear 
behaving  as  a  random  rough  surface.  In  the  linear  approximation,  the  statistics  of 
the  ocean  surface  are  homogeneous,  stationary,  and  Gaussian.  With  these  assump¬ 
tions,  several  fundamental  electromagnetic  backscatter  models,  in  good  agreement 
with  observations,  have  been  derived  (Barrick  and  Peake,  1968;  Barrick,  1968,  1972; 
Ilasselmann,  1971;  Kirshen,  1971;  Daley,  1973;  Valenzuela,  1974,  1978;  among  oth¬ 
ers).  These  models  are  based  upon: 

•  The  high  frequency  methods,  including  geometric  and  physical  optics  which 
yield  specular  point  models  for  smooth  and  undulating  surface  (Kerr,  1951). 

•  The  low  frequency  method,  using  the  perturbation  method  for  a  slightly  rough 
surface  yields  the  Bragg  scattering  model  (Rice,  1951). 

•  Combination  of  high  and  low  frequency  methods  renders  the  two-scale  com¬ 
posite  model  for  rough  surfaces. 

The  purpose  of  these  fundamental  models  of  microwave  scattering  from  the 
ocean  surface  is  to  convert  the  received  electric  signal  at  the  sensor  into  a  mean¬ 
ingful  geophysical  signal.  These  models  (specular,  Bragg-Rice,  and  composite)  are 
described  extensively  in  a  review  paper  by  Valenzuela  (1978),  which  remain  the  ba¬ 
sis  for  studies  of  ocean  surface  dynamics  sensed  by  microwave  sensors  even  in  the 
present  time  (see  Figure  2.1).  Recently,  the  composite  model  has  been  extended  to 
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incorporate  a  third  intermediate  scale  utilized  specifically  for  SAR  imagery  appli¬ 
cations.  This  theory  is  summarized  in  great  detail  in  a  consensus  paper  by  Has- 
selmann,  Raney,  Plant,  Alpers,  Shuchmann,  Lyzenga,  Rufenach,  and  Tucker,  1985. 
The  scattering  of  electromagnetic  fields  from  random  rough  surface  differs  from  de¬ 
terministic  surfaces  scattering  by  being  statistical  in  nature.  Mean  values  of  the 
fields  and  power  may  be  determined  by  ensemble  averaging  these  quantities  over 
the  assumed  or  known  statistics  of  the  rough  surface  (Beckmann  and  Spizzichino, 
1963). 


Incident 

Wave 


Slightly  rough 


Incident  Incident 

Wave  Wave 


Undulating  surface  Rough  surface 


Figure  2.1:  Illustration  of  EM  scattering  from  ocean  surface. 

Following  Valenzuela  (1978),  the  scattered  power  density  is  given  by  the  mean 
magnitude  of  the  Poynting  vector,  which  is  a  cross  product  of  the  electric  and 
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magnetic  fields 


P,  =  l-R,(E,  X  h;)  (2.1) 

where  H*  is  the  complex  conjugate  of  the  scattered  magnetic  field  of  the  EM  wave, 
and  E4  is  the  scattered  electric  field  of  the  EM  wave.  The  brackets  denote  ensemble 
average. 

A  measure  of  scattering  power  from  a  rough  surface  can  be  related  to  the 
radar  cross-section,  which  is  proportional  to  the  ratio  of  the  scattered  power  from 
the  surface  to  the  incident  power.  For  the  far  field: 

<r  =  Urn  (toR'm’l/m2)  (2.2) 

where  |£i|2  and  {|E4|2)  are  the  incident  and  mean  scattered  power  respectively.  R 
is  the  distance  of  the  observation  point  from  the  scattering  rough  surface. 

In  the  following  sections  the  specular,  Bragg-Rice,  and  the  composite  scat¬ 
tering  models  are  presented. 


to  radar 


Figure  2.2:  Specular  point  backscatter  (from  Stewart,  1985). 

1.  Specular  Point/Physical  Optics  Model 

In  this  model,  the  reflections  are  from  surface  slopes  which  are  oriented 
normal  to  the  line  of  sight  of  the  radar  (Hagfors,  1964).  The  principal  assumptions 


8 


of  this  model  are  that  the  surface  radius  of  curvature  is  almost  everywhere  much 
larger  than  the  radar  wavelength  (physical  optics  or  Kirchhoff  method),  and  that 
multiscattering  and/or  shadowing  are  negligible  (see  Figure  2.2).  The  range  of 
validity  of  physical  optics,  or  the  tangent  plan  approximation,  was  shown  by  Shmelev 
(1972)  as: 

2 rkji  cos  7  >  1  and  - rkR  cos3  7  »  1  (2.3) 

where  r  is  the  local  radius  of  curvature  of  the  surface,  7  and  kn  are  the  incidence 
angle  and  the  wavenumber  of  the  electromagnetic  wave  in  free  space  respectively. 

The  backscatter  power  from  a  single  facet  (local  tangent  plan)  is  a  func¬ 
tion  of  the  probability  that  the  tangent  plane/facet  will  be  oriented  perpendicularly 
to  the  incidence  direction  of  the  electromagnetic  wave.  Summing  up  the  contribu¬ 
tions  from  all  properly  oriented  facets  within  the  resolution  cell  of  the  radar  results 
in  the  total  backscatter  specular  power. 

Applying  the  physical  optics  approximation,  Barrick  (1968)  derived  a 
specular  model  to  predict  the  normalized  radar  cross  section  for  backscatter  from  a 
finitely  conductive  rough  surface: 


<70  =  (?r/ cos' 


(2.4) 


where  p{px,r]y)  is  the  joint-probability  density  of  orthogonal  slopes  tjx  and  r)v  of  the 
sea  surface  rj  with  respect  to  the  incident  EM  wave  and  p(0)  is  the  Fresnel  reflection 
coefficient  for  normal  incidence.  The  slopes  are  the  horizontal  derivatives  of  the 
local  vertical  displacements  of  water  surface  from  mean  level. 

A  simplified  form  of  Equation  (2.4)  is  achieved  by  assuming  that  the 
ocean  surface  roughness  has  a  Gaussian  random  distribution.  Thus,  for  an  upwind 
looking  radar,  Equation  (2.4)  becomes: 

a°  =  2S^5~ sec^0)  exp  ta°2  °/2SD  (2.5) 
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where  Si  and  S%  are  the  mean  square  slopes  of  upwind  and  crosswind  respectively. 
Verified  by  observations,  Equation  (2.5)  seems  to  be  a  good  approximation  for  near 
normal  backscatter  from  an  undulating  surface. 

Radar  observations  for  the  ocean  indicate  that  the  specular  backscatter 
is  2  to  3  times  larger  from  a  trough  than  from  a  crest  (Yaplee,  1971;  Valenzuela, 
1978).  This  could  be  due  to  “compressibility”  of  the  short  gravity-capillary  waves 
at  the  crests  of  the  long  waves,  resulting  in  a  diffuse  scatter  which  reduces  the 
specular  backscatter.  The  reflection  properties  of  a  surface  are  known  to  be  reduced 
by  increase  in  roughness.  Away  from  normal  incidence,'  the  backscatter  is  no  longer 
specular,  but  diffuse,  being  produced  by  shorter  waves  which  have  steeper  slopes 
(roughness). 

2.  The  Bragg-Rice  Model 

For  a  slightly  rough  surface,  to  first  order,  the  backscatter  energy,  in 
directions  of  at  least  20  degrees  from  normal  incidence,  is  proportional  to  the  ripple 
•variance  spectrum  evaluated  at  the  Bragg  resonant  wave  number  (Hasselmann  et 
al.,  1985) 

fc  =  2fcflsin7  (2.6) 

where  k  is  the  ocean  wavenumber,  and  kn  sin  7  is  the  projection  of  the  incident  radar 
wave  number  onto  the  ocean  surface.  This  diffraction  grating  (resonant)  backscat- 
tering  theory  for  e.m.  waves,  was  first  developed  by  Rice  (1951)  and  demonstrated 
experimentally  by  Crombie  (1954)  with  an  IIF  Doppler  radar.  The  Doppler  spec¬ 
trum  showed  two  peaks  shifted  from  the  EM  frequency  by  the  positive  and  negative 
Bragg  waves  frequency,  for  waves  propagating  towards  and  away  from  the  sensor. 

Since  the  resonant  (Bragg-Rice)  backscattering  mechanism  has  been 
shown  to  be  the  dominant  backscatter  for  the  ocean  surface  for  incidence  angles 
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ranging  from  20  to  70  degrees,  typical  for  SAR  remote  sensing,  this  backscattering 
mechanism  is  reviewed  in  more  detail,  together  with  the  composite  surface  model 
introduced  in  the  next  section. 

The  main  assumptions  in  deriving  the  first  order  Bragg- Rice  model  are: 

•  Local  matching  of  the  fields  at  the  surface  with  linearized  boundary  conditions. 
Hence,  the  surface  vertical  perturbations,  A,  of  the  Bragg  ripples  are  small  with 
respect  to  the  radar  wavelength  (Rayleigh  condition). 

•  The  number  of  resonant  scatterers  is  large  enough  to  allow  a  proper  statistical 
description  of  the  backscattered  reflectivity. 

•  The  slopes  of  the  Bragg  waves  are  assumed  small.  Jordan  and  Lang  (1972) 
argue  that  this  condition  represents  the  real  ocean  because  water  waves  break 
long  before  they  become  large  enough  to  violate  the  linear  boundary  condition. 


Based  on  these  assumptions,  the  resonant  scatter  theory  yields  a  nor¬ 
malized  radar  cross  section  (Wright,  1966;  Barrick,  1972;  Jordan  and  Lang,  1979): 

<70  =  87r^|<7O(7)|V(2fc/{sin7,0)  (2.7) 


where  0,7(7)  are  the  first  order  scattering  coefficients  and  t/>(&x,  ky )  is  the  two  dimen¬ 
sional  ocean  wave  number  spectrum  properly  normalized.  The  indices,  (i,j)  denote 
transmitted  and  received  polarizations  respectively: 


9ij{  7)  = 


9vh{i) 


9hv{  7) 
9vv{l) 


(2.8) 


Depolarization  effects  for  the  ocean  surface  are  generally  smaller  (of  higher  order) 
than  the  copolarized  contribution,  therefore  only  scattering  coefficients  with  un¬ 
changed  polarization  with  respect  to  the  transmitted  polarization  are  considered 
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here: 


„  /  X  _  (£r  ~  1)  COS2  7 

9HH(V  —  r  .  1  1  2 

cos  7  +  (£r  -  sin  7}sJ 

.  .  (sr  -  l)[er(l  +  sin2  7)  -  sin2  7]  cos2  7 

gwiV  -  r  .  ,  xll2 

ler  cos  7  +  (er  -  sm  7}  *  j 

where  er  is  the  relative  complex  dielectric  constant. 


gwi'r)  = 


(2.10) 


Away  from  normal  incidence,  Bragg  resonant  scattering  yields  an  in¬ 
creasing  polarization  ratio  with  increasing  incidence  angle,  with  a  larger  vertical  po¬ 
larization  scattering  coefficient.  For  a  surface  with  perfect  conductivity,  er  =  —too, 
which  is  a  reasonable  assumption  for  the  ocean  surface  for  lower  radar  frequencies 
(i.e.,  IIF),  Equations  (2.9)  and  (2.10)  simplify  (Stewart,  1985): 


9uh{0)  ~  cos2  6 


(2.11) 


gvv{0)  ~  1  +  sin2  0  (2.12) 

Therefore,  for  increasing  incidence  angles,  the  horizontal  polarization  Bragg  scat¬ 
tering  coefficient  decreases  in  magnitude,  and  is  more  sensitive  to  incidence  angles 
relative  to  vertical  polarization  scattering  coefficient  except  for  near  grazing  inci¬ 
dence  where  the  trapped  waves  and  shadowing  become  dominant. 

The  backscatter  resonant  cross-section  is  directly  proportional  to  the 
two-dimensional  ocean  wave  number  spectrum  satisfying  the  Bragg  resonant  condi¬ 
tion  (Equation  2.7),  which  in  turn  suggests  the  usefulness  of  microwave  backscatter 
for  probing  ocean  waves.  If  the  Bragg  water  waves  are  not  the  dominant  waves  of  the 
ocean,  higher  order  Doppler  shifted  backscatter  will  occur  whenever  the  summation 
of  wavenumber  pairs  of  long  waves  meet  the  Bragg  resonant  condition  (Hasselmann, 
1966): 

2&/*sin  6  =  ±ki  -f  (±^2)  (2.13) 
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a>£)  =  ±u>i  -f  (dhu>2)  (2.14) 

where  k\  and  &2  are  the  wave  numbers  of  the  waves.  According  to  perturbation 
mathematical  constraints,  the  second  order  contribution  should  be  smaller  compared 
with  the  first  order  Bragg  contribution.  In  practice,  some  of  these  higher  order 
contributions  representing  initial  manifestation  of  tilting  and  multiple  scattering 
can  become  comparable  to  the  first  order  contribution  or  larger.  Typical  relative 
amplitude  Doppler  spectra  were  observed  by  Barrick  et  al.  (1974)  using  a  HF 
frequency  Doppler  radar  (see  Figure  2.3).  The  Doppler  spectrum  in  Figure  2.3  also 
shows  a  frequency  shift  of  the  first  order  Bragg  lines  due  to  advection  of  the  waves 
by  currents. 


Figure  2.3:  Bragg  phenomena  in  Doppler  frequency  domain  (from  Valen¬ 
zuela,  1978). 

Wind  may  induce  surface  advection  drift,  which  will  add  to  the  Doppler 
shift  of  the  phase  speed.  Also,  rotational  (finite  amplitude)  effects  may  modify 
the  phase  speed  of  short  waves.  At  high  winds  as  the  ocean  becomes  rougher,  the 
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higher  order  contributions  become  significant,  and  the  Doppler  spectra  becomes 
broader.  The  broader  Doppler  spectra  occur  when  the  wind  dominant  waves  are 
longer  than  the  Bragg  waves  on  the  ocean.  When  the  Bragg  resonant  waves  are 
generated  on  longer  dominant  ocean  waves,  the  mean  surface  for  the  Bragg  waves 
is  no  longer  flat  but  undulating.  For  this  typical  ocean  surface,  a  composite,  or 
two-scale  hydrodynamic  model  is  suggested  in  the  literature. 

3.  Composite  Surface  Model 

The  observation  that  Bragg  resonant  theory  by  itself  does  not  predict  in¬ 
creased  radar  backscatter  in  the  presence  of  longer  waves  led  Bass  et  al.  (1968)  and 
Wright  (1976)  to  suggest  a  more  tractable  two  scale-composite  surface  backscatter 
model.  This  model  represents  short  Bragg  waves  riding  on  longer  waves  by  ap¬ 
proximating  the  ocean  surface  as  a  sum  of  infinite  number  of  slightly  rough  patches 
(Figure  2.4).  The  composite  model  predicts,  therefore,  the  radar  cross-section  for  a 
single  patch  as  given  by  (Wright,  1968;  Valenzuela,  1978;  Brown,  1978): 

°Hu{0i)  =  8xk4n  (j‘C°S^')  ghh(9i)  +  -j  aM)  il>(2ka,2kj  sin  tf)  (2.15) 

where  a;  =  sin  7,-,  7 ,■  =  cos(7  -f  ^)cos6,  a  =  sin(7  +  ty)  and  7 ;  =  cos  0,-,  and 
are  the  local  sea  surface  slopes  (small  angles),  ghh{0i)  and  gvv{0i)  are  the  first  order 
scattering  coefficients  of  a  local  patch,  oyv  is  computed  with  a  similar  formula,  but 
it  is  much  less  sensitive  to  variations  in  local  facet  tilts.  Thus,  the  contribution  from 
all  slightly  rough  patches  over  the  distribution  of  slopes  p{r}x,Vv)  °f  the  dominant 
waves  within  the  radar  resolution  cell,  is  given  by: 

/OO  f  OO 

/  07///(#,)p(tan$,tan<5)<i(tan'I')d(tan<5)  (2.16) 

•OO  J  —  OO 

Normally,  the  backscatter  from  a  slightly  rough  surface  modulated  by 
long  waves  is  approximated  including  both  the  contribution  from  a  specular  point 
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lo  radar 


Figure  2.4:  Illustration  of  composite  surface  approximation  (from  Stew¬ 
art,  1985). 

model  and  that  predicted  by  the  composite  surface  model  (Hasselmann  et  al,  1985). 
The  tilted  facets  significantly  modify  the  radar  cross-section,  especially  for  the  hor¬ 
izontally  polarized  radiation  with  a  depolarized  contribution  (Stewart,  1985). 

Towards  grazing  incidence  the  polarization  ratio  decreases  as  the  sea  ap¬ 
pears  rougher,  a  phenomenon  which  illuminated  Wright  to  develop  the  composite 
model.  The  cross-section  polarization  ratio  as  a  function  of  incidence  angle  and 
roughness  parameter  is  illustrated  in  Figure  2.5.  At  grazing  incidence  a  modified 
model  including  a  combination  of  composite  surface  and  wedge  scattering  was  sug¬ 
gested  by  Kalmykov  and  Pustovoytenko  (197C)  to  explain  Doppler  spectra  from  the 
ocean.  This  combined  backscalter  model  predicts  larger  backscatter  and  greater 
Doppler  shift  for  horizontal  polarization  (at  grazing  incidence),  rather  than  for  ver¬ 
tical  polarization. 

4.  Application  of  Composite  Model  to  SAR 

For  purposes  of  interpreting  coherent  imaging  observations,  the  actual 
scattering  from  a  moving  surface  like  the  ocean  can  be  represented  by  a  superposi- 
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Figure  2.5:  Cross-section  ratio  due  to  polarization  (from  Valenzuela, 
1978). 


tion  of  the  specular  and  composite  theories  as  described  above  (however,  it  is  not 
uncommon  to  include  specular  scattering  in  the  definition  of  composite  theories), 
plus  complex  motion  effects  from  a  non-rigid  distributed  surface.  The  interpretation 
of  images  from  the  ocean  surface,  is  still  an  area  of  active  research.  Elachi  (1977); 
Hasselmann  et  al.  (1985),  and  Alpers  and  Bruening  (1986)  have  given  some  insight 
to  SAR  imaging  models  of  the  ocean  surface.  Some  of  their  considerations  will  be 
described  in  Chapter  III.  D. 
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III.  PRINCIPAL  MECHANISMS  OF  SAR 
AND  INSAR  IMAGERY  RELEVANT  TO  THE 

OCEAN  SURFACE 

A.  INTRODUCTION 

In  this  chapter  the  fundamental  principles  of  SAR  imaging  are  briefly  de¬ 
scribed  first,  since  the  suggested  interferometric  technique  is  based  on  conventional 
SAR.  A  synthetic  aperture  radar  is  an  imaging  sensor,  with  a  high  (two  dimensional) 
spatial  resolution  in  azimuth  (along  track)  and  in  range  (across  track)  directions 
(Ilarger,  1970;  Ilovanessian,  1979).  High  azimuthal  resolution  is  achieved  by  trans¬ 
lating  a  coherent  radar  on-board  a  moving  platform,  storing  the  amplitude  and  phase 
histories  returned  from  the  target,  and  subsequently  processing  the  stored  data  to 
remove  known  quadratic  phase  errors  across  the  antenna  beam.  Improved  range 
resolution  is  achieved  by  using  spread  spectra  techniques,  such  as  pulse  compression 
or  bi-phase  modulation  of  the  radar  carrier  frequency  (Skolnik,  1980). 

The  range  scanning  is  at  the  speed  of  light  in  air,  but  the  azimuthal  scanning 
is  at  the  speed  of  the  carrier  platform  (Ilasselmann  et  al.,  1985).  For  processing 
purposes  this  fact  leads  to  substantial  decoupling  of  the  range  and  azimuth  process¬ 
ing.  The  azimuth  channel  is  related  directly  to  the  dynamic  distortions  associated 
with  imaging  of  ocean  surface  features,  because  the  time  scales  of  the  system  and 
the  imaged  wave  fields  are  of  the  same  order.  Inherently,  the  range  channel  does  not 
have  such  problems.  In  this  chapter,  some  fundamental  parameters  characterizing 
SAR  and  INSAR  operation  are  described. 
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B.  FUNDAMENTAL  PARAMETERS 


1.  Azimuth  Resolution 


At  any  given  time  and  range,  the  contributions  to  the  radar  image  from 


each  stationary  scattering  point  at  that  range  can  be  identified  uniquely  by  its 


Doppler  shift,  caused  by  the  sensor  motion.  This  allows  the  summation  of  different 


contributions  from  different  scattering  points  to  produce  a  phase  history  in  the 


aperture  synthesis. 


An  antenna  with  a  real  aperture  D  has  an  azimuthal  coverage  of  (see 


Figure  3.1) 


ACROSS 

THE 

TRACK 


ILLUMINATED 
PATCH  OP 
SURFACE 


Figure  3.1:  SAR  geometry. 


(3.1) 


where  R  is  the  radar  target  slant  range,  and  A  is  the  radar  wavelength.  A Rx  is  the 


real  azimuth  resolution  with  antenna  of  size  D.  The  length  of  aperture  which  can 


be  synthesized  depends  on  the  distance  the  platform  travels  while  a  point  scatterer 


is  still  illuminated.  The  synthetic  aperture  for  the  geometry  in  Figure  3.1  is  A Rx, 
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and  its  maximum  azimuth  resolution  for  a  stationary  rigid  point  target  is: 


A  R  D  .  . 

pazG  =  2A Tx  =  2  (3-2) 

the  factor  of  2  enters  because  the  effective  synthetic  aperture  “beam-width”  is  half 
of  a  real  aperture.  From  Equation  (3.2),  the  synthetic  aperture  resolution  is  inde¬ 
pendent  of  the  wavelength  and  range  of  the  radar,  and  increases  with  a  decrease 
in  the  physical  size  of  the  antenna  (Hovanessian,  1979;  Robinson,  1984;  Stewart, 
1985). 

2.  Phase  History  (of  stationary  point  target) 

Since  the  resolution  in  azimuth  direction  is  based  on  the  local  Doppler 
shift  due  to  aircraft  motion,  which  in  turn  is  the  rate  of  change  of  phase,  the  phase 
history  is  defined,  which  is  the  principle  of  synthetic  aperture  imaging  radar.  Con¬ 
sider  a  moving  SAR  with  velocity  V  as  illustrated  in  Figure  3.2.  By  the  time  t,  the 
radar  moved  a  distance  L=Vt,  and  the  range  to  a  stationary  point  target  changed 
from  R\  to  Ro  to  Rx.  Since  the  phase  of  the  received  signal  is  proportional  to  the 
distance  traveled,  the  phase  of  this  signal  will  vary  with  time.  The  range  R  can  be 
written  as: 

R  =  [1%  +  (Vt)2}1'2  (3-3) 


where  Rq  is  target  range  of  closest  approach.  By  a  Maclaurin  series  expansion 
around  t  =  0  Equation  (3.3)  can  be  rewritten  as: 


R=^k{Vt?-i kVi)'+- 


(3.4) 


Neglecting  higher  order  terms,  the  phase  shift  of  a  stationary  point  target  is  (Hov¬ 
anessian,  1979): 


A$(<)  = 


2tt  {Vi)2 
A  Ro 


(3.5) 
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Figure  3.2:  Top  view  of  range  to  a  point  target  as  a  function  of  aircraft 
velocity  (from  Hovanessian,  1979). 


Equation  (3.5)  states  that  the  Doppler  shift  from  a  stationary  target  is  a  linear 
function  of  time  (chirp).  Two  additional  fundamental  parameters  which  determine 
SAR  performance  are  the  processor  coherent  integration  time  T  and  the  Doppler 
bandwidth  D. 

3.  Coherent  Integration  Time 

The  upper  limit  of  the  theoretical  coherent  integration  time  is  deter¬ 
mined  by  the  real  aperture  azimuth  resolution  or  more  precisely,  by  the  real  antenna 
beamwidth: 


T  = 


&RX  _  PR 
V  '  V 


(3.6) 


where  ft  is  the  real  aperture  beamwidth. 

4.  Doppler  Bandwidth 

Consider  the  scenario  depicted  in  Figure  (3.3).  The  Doppler  shift  of  a 
stationary  point  target  B  is: 


Ib  - 


2V  cos  0 

A 


(3.7) 
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Figure  3.4:  Azimuth  image  shift. 


C.  EFFECT  OF  POINT  TARGET  MOTION 

Raney  (1971)  and  Alpers  and  Bruening  (1986),  among  others,  give  a  clear 
description  of  the  effect  of  point  target  motion  on  SAR  image.  Following  Alpers 
and  Bruening  (1986),  consider  a  point  target  with  radial  velocity  Ut  in  the  cross 
track  (range)  direction  as  depicted  in  Figure  3.5.  The  SAR  carrier  is  moving  at  V 
velocity  and  R  is  the  distance  between  the  SAR  antenna  and  the  point  target. 

A  SAR  which  is  tuned  for  stationary  targets,  locates  targets  in  azimuth  at 
zero  Doppler  frequency.  Hence,  a  stationary  reasonable  target  will  be  located  at 
right  angles  to  the  carrier  path.  However,  for  the  moving  target  (Figure  3.4),  the 
zero  Doppler  frequency  position,  will  be  shifted  in  azimuth  at  an  angle  a  = 
relative  to  a  position  of  a  stationary  target.  Thus  the  azimuth  image  shift  of  a 
target  with  radial  velocity  is: 

A  x  =  -^Ur  (3.10) 

Evidence  of  this  shift  is  found  in  SAR  imagery  for  moving  rigid  targets  with  velocity 
component  in  the  range  direction,  like  ships  or  trains.  Their  images  are  shifted  in 
azimuth  direction  compared  to  their  wakes  or  tracks.  Keeping  in  mind  fundamental 
principles  of  SAR,  the  imaging  of  the  moving  ocean  surface  is  discussed  next. 


Figure  3.5:  Radar  hydrodynamic  and  straining  modulation  transfer  func¬ 
tion  and  SAR  velocity  bunching  azimuthal  transfer  function  as  a  function 
of  azimuth  angle  (from  Alpers  et  al.,  1981). 

D.  APPLICATIONS  TO  IMAGING  OF  MOVING 
OCEAN  SURFACE 

As  described  in  section  B  of  this  chapter,  microwave  scattering  from  the  ocean 
surface  for  conventional  radars  is  modeled  to  an  acceptable  extent.  However,  when 
dealing  with  imaging  radars,  two  new  parameters  improving  the  radar  performance 
and  in  a  way  complicating  the  backscattering  process  are  introduced: 


t  High  spatial  resolution  (order  of  tens  of  meters). 
•  Doppler  information  for  image  generation. 
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The  full  implication  of  these  two  additional  factors  on  the  generation  of  ocean  sur¬ 
face  imagery  is  still  in  active  research  (Elachi  and  Brown,  1977;  Kasilingham  and 
Shemdin,  1988). 

For  a  stationary  surface,  a  SAR  image  is  a  spatial  two-dimensional  presenta¬ 
tion  of  the  variability  of  the  local  coherent  backscatter  cross  section  of  the  surface 
(Elachi  and  Brown,  1977).  For  the  ocean  surface  which  is  nonstationary,  SAR  image 
is  a  map  of  the  covariance  of  complex  reflectivity  (Hasselmann  et  al.,  1985;  Kasil¬ 
ingham  and  Shemdin,  1988).  Due  to  the  difficulties  in  interpreting  SAR  imagery 
from  the  ocean  surface,  two  major  resolvable  long  wave  imaging  mechanism  are 
still  recognized,  though  they  represent  conflicting  view  points.  One  theory  is  the 
velocity  bunching  approach,  which  is  based  on  the  coherent  particle  orbital  velocity 
of  the  dominant  wavefield  (Elachi  and  Brown,  1977;  Alpers  and  Rufenach,  1979; 
Swift  and  Wilson,  1979;  Keller  and  Wright,  1983;  Hasselmann  et  al.,  1985;  Alpers 
and  Bruening,  1986).  The  other  theory  is  the  distributed  surface  theory  which  pre¬ 
dicts  that  the  SAR  image  is  composed  of  amplitude  modulation  due  to  straining  of 
short  backscaltering  waves  and  a  phase  modulation  due  to  tilting  by  long  waves. 
It  assumes  that  the  ocean  surface  is  correlated  over  distances  larger  than  the  SAR 
intrinsic  resolution.  (Jain,  1978;  Harger,  1986). 

1.  Velocity  Bunching  Imaging  Mechanism 

The  velocity  bunching  theory  is  based  on  the  fundamental  assumption 
that  the  scattering  surface  can  be  divided  into  individual  and  uncorrelated  scattering 
elements  whose  dimensions  are  small  compared  to  the  SAR  resolution  cell  (~  10—20 
m)  and  large  compared  to  the  radar  wavelength  (5-7  times).  Each  scattering  element 
(facet)  is  associated  with  the  local  orbital  velocity  and  acceleration  of  the  long 
underlying  waves.  This  facet  motion  gives  rise  to  alternating  concentrating  and 
spreading  of  the  position  of  the  scattering  facets  causing  amplitude  modulation 
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phase  locked  to  the  dominant  waves.  Following  Vachon  et  al.  (1988),  the  radial 
line  of  sight  velocities  of  a  monochromatic  long  wave  associated  with  a  scattering 
element  is  given  by 

U(x,  t)  =  w/l(sin27sin2  <f>  +  cos27)1/2sin(k  •  x  —  wt  -f  (£<77  cos  <f> ))  (3.11) 

where  <f>  is  the  angle  between  the  aircraft  velocity  vector  and  the  wavenumber  vector 
k  on  the  ocean  surface  plane,  and  A  is  the  wave  amplitude.  For  swell  with  a 
period  much  longer  than  the  SAR  integration  time,  a  scattering  element  originally 
at  azimuthal  coordinate  x  will  be  shifted  to  (Raney,  1988) 

x'  —  x  +  c'  sin(fccos  <f)x)  (3.12) 

A  constant  phase  term  was  not  included  in  Equation  (3.12).  c'  is  the 
fractional  amplitude  of  the  azimuthal  image  shift  of  a  scattering  element  and  governs 
the  mapping  character,  c'  is  defined  as  (Alpers,  1983) 

c'  =  (L2)  ^/i(sin2  7sin2  <f>  +  cos2  7)1/2  (3.13) 

Under  certain  situations  the  processed  SAR  image  can  appear  with 
bright  intensity  patterns  which  have  twice  the  spatial  period  of  the  actual  long 
wave  on  the  ocean,  which  can  confuse  the  interpretation  of  ocean  wave  analysis. 
This  critical  situat  ion  is  quantified  by  the  parameter  c  and  becomes  critical  for 

c'r  =  L/2ncos<f>  (3-14) 

for  c'  >  L/2ncos<j)  the  scattering  element  redistribution  becomes  many-to-one  and 
results  in  two  bright  high  reflectivity  bands  associated  with  one  spatial  period  of 
the  actual  ocean  wave.  For  c'  <  /y/27rcos<^,  the  scattering  element  redistribution 
is  one  to  one  and  results  in  a  single  unambiguous  high  reflectivity  band  related  to 
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each  spatial  period  of  the  actual  ocean  wave.  The  parameter  c'r  is  defined  by  Alpers 
and  Rufenach  (1979)  as  a  measure  of  transition  between  linear  and  nonlinear  SAR 
imaging  (see  Figure  3.4).  In  terms  of  the  significant  ocean  wave  height,  the  critical 
imaging  occurs  when 


/M3/2 _ V _ 

\  jr  /  g1/2R  cos  <^(sin2  7  sin2  <f>  +  cos2  7)1/2 


(3.15) 


In  the  above  velocity  bunching  model  the  ocean  surface  wets  assumed  commonly 
(Vachon  et  al.,  1988)  to  be  invariant  in  its  structure  and  scatter  dynamics  during 
the  SAR  image  formation  period  (coherent  integration  time).  The  integration  time  is 
considered  to  be  much  less  than  the  dominant  imaged,  wave  period.  This  might  not 
be  completely  true,  especially  for  airborne  SAR.  Kasilingam  and  Shemdin  (1988) 
and  others  suggest  that  the  noncoherent  scene  motion  (phase  velocity)  over  the 
image  formation  period  should  be  explicitly  considered  in  the  velocity  bunching 
model. 


Additional  motion  effects  of  the  scattering  elements  that  are  sensed  co¬ 
herently  by  the  SAR  radar  is  the  dominant  wave  acceleration  that  corresponds  to 
the  orbital  motion  and  the  scattering  element  decorrelation  time  scale.  Together 
with  the  translation  of  wave  patterns  during  the  integration  time  which  enters  the 
image  noncoherently,  these  additional  ocean  surface  dynamics  blurs  the  image.  The 
scene  coherence  time  effect  is  discussed  in  detail  in  Chapter  IV.  The  acceleration 
defocus  is  discussed  in  Chapter  VIII. 

2.  Noncoherent  Motion  Effects 

The  translation  of  mean  wave  patterns  of  reflectivity  during  the  course 
of  integration  time  enters  noncoherently  into  the  image  in  two  ways:  through  distor¬ 
tion  in  wave  length  and  direction  known  as  scanning  distortion  and  through  image 
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smearing  (defocus),  also  known  as  look  misregistration  for  a  multilook  imaging  sys¬ 
tem. 

For  an  monochromatic  imaged  swell,  the  noncoherent  smearing  can  be 
corrected  by  adjusting  the  processor  focus  perturbation  to  one  half  the  projected 
wave  phase  velocity.  This  operation  will  enhance  the  SAR  ocean  wave  image  (Jain 
and  Shemdin,  1983). 

Scanning  distortion  is  associated  with  imaging  a  moving  scene  from  a 
moving  sensor.  The  along  track  motion  of  the  scattering  cross-section  pattern  mov¬ 
ing  at  the  long  wave  phase  velocity  distorts  the  azimuthal  imaged  wavenumber. 
This  noncoherent  distortion  is  not  related  to  the  coherent  phase  history  of  individ¬ 
ual  scattering  elements.  If  the  dominant  imaged  wave  is  propagating  only  in  the 
cross-track  direction,  only  a  distortion  in  wave  direction  will  be  practically  mani¬ 
fested.  The  degree  of  distortion  in  either  wave  orientation  is  proportional  to  the 
ratio  between  the  wave  phase  velocity  and  platform  velocity  C/V\  The  distorted 
azimuthal  wavenumber  is  taken  to  be 

k'az  =  kaz±k^  (3.16) 

The  estimation  of  scanning  distortion  from  two  opposite  flights  provides 
an  opportunity  to  resolve  the  ambiguity  of  180°  in  propagation  direction  that  is 
inherent  to  remote  sensors. 

3.  SAR  Two-Scale  Backscattering  Model 

Ilasselmann  et  al.  (19851  gives  a  rather  complete  theory  for  Bragg  scat¬ 
tering  using  the  two  scale  approach  in  a  review  paper.  They  distinguish  between 
two  types  of  two-scale  models: 

•  The  electromagnetic- hydrodynamic  two- scale  model. 
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•  SAR  two-scale  model. 


The  basic  difference  between  these  models  is  the  scale  of  separation  between  the  two 
wave  regimes.  In  both  models,  the  wave  field  is  treated  statistically  for  wave  numbers 
larger  than  the  scale  of  separation  and  deterministically  for  lower  wave  numbers.  In 
the  electromagnetic-hydrodynamic  two-scale  model,  the  scale  of  separation,  fcem/,,  is 
defined  by  the  expansion  requirements  of  electromagnetic  hydrodynamic  theory.  The 
scale  of  separation  has  to  be  chosen  one  order  of  magnitude  smaller  than  the  Bragg 
backscattering  wavenumber,  h.  On  the  other  hand  in  the  SAR  two-scale  model, 
the  scale  of  separation,  k>ar ,  is  determined  by  the  inverse  of  SAR  resolution  scale. 
The  SAR  two-scale  imaging  model  is  based  on  the  physics  of  the  more  fundamental 
two-scale  electromagnetic-hydrodynamic  model.  Typically  the  fcem/,  is  one  order  of 
magnitude  larger  than  k3ar  (see  Chapter  IV). 

4.  SAR  Mapping  Transfer  Function 

In  addition  to  the  Doppler  induced  modulations  (velocity  bunching)  typ¬ 
ical  for  coherent  imaging  radars,  real  radar  cross-section  modulation  tilt  and  hydro- 
dynamic  modulation  are  commonly  used  to  describe  ocean  surface  imaging.  These 
imaging  mechanisms  govern  the  imaging  process  for  waves  propagating  close  to  range 
(cross-track)  direction.  In  this  orientation  the  velocity  bunching  is  nil  (see  Figure 
3.4).  For  moderate  sea  states  the  mapping  of  the  wavefield  might  be  described  by  a 
linear  transfer  function 

Rsar  =  +  RhydT 

Usually  this  situation  does  not  occur,  since  the  mapping  is  non-linear  and  the  super¬ 
position  principle  does  not  hold,  i.e.,  the  image  spectrum  is  distorted  from  the  actual 
wavefield.  Electromagnetic  interaction ,  or  tilt  modulation  of  backscatter  radar  cross 
section,  is  caused  by  the  actual  variation  of  local  incidence  angle.  Hydrodynamic 


_j_  ^bunching 


(3.17) 
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modulation  of  backscatter  radar  cross  section  is  caused  by  the  interaction  of  the 
Bragg  waves  with  the  long  waves  on  which  they  are  superimposed.  The  long  wave 
orbital  velocity  strains  the  surface  producing  areas  of  current  convergence  and  di¬ 
vergence,  resulting  in  hydrodynamic  modulation  of  the  amplitude  and  wavenumber 
of  the  Bragg  waves.  The  surface  displacement  distorts  the  air  flow  over  the  sur¬ 
face,  causing  short  waves  to  grow  at  varying  rates  (modulated  wind  stress)  for  a 
given  wind  speed,  depending  on  their  relative  position  on  the  long  wave  (private 
communication,  Valenzuela,  1989). 

E.  INTERFEROMETRIC  SAR  (INSAR) 

1.  Principle  of  Imaging  of  Ocean  Surface  Velocity  Fields 

Interferometric  synthetic  aperture  radar  (INSAR)  is  a  relatively  new 
method  for  measuring  ocean  surface  wave  fields,  spectra,  and  currents.  Unlike 
conventional  SAR,  INSAR  has  the  potential  for  measuring  unambiguously  virtual 
surface  velocity  fields.  INSAR  combines  large  coverage  with  high  nominal  spatial 
resolution  (~  12  m  in  this  work),  good  sensitivity  and  appears  suitable  for  both 
airborne  and  probably  satellite  implementation.  The  airborne  INSAR  technique 
has  been  developed  at  Jet  Propulsion  Laboratories  (JPL)  by  Goldstein  and  Zebker 
(1987),  and  was  installed  first  on  a  NASA  CV990  aircraft  and  recently  on  a  DC-8 
for  evaluation  experiments. 

Following  Goldstein  and  Zebker  (1987),  the  basic  technique  consists  of 
aligning  two  antennas  along  the  radar  flight  path.  Pulses  are  transmitted  from  the 
rear  antenna  and  are  received  simultaneously  by  both  antennas.  The  resulting  sig¬ 
nals  are  converted,  separately,  into  two  complex  images  in  the  conventional  manner 
for  SAR.  Any  line  of  sight  motion  of  the  identical  imaged  surface  (resolution  cell) 
observed  by  both  antennas  is  translated  into  phase  difference  when  these  resolution 
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cells  are  combined  interferometrically.  This  phase  difference,  ct,  is  proportional  to 
the  virtual  measured  radian  Doppler  shift,  A u>£>,  that  occurs  in  the  time  it  takes  for 
the  rear  antenna  to  move  to  a  place  of  symmetry  with  respect  to  the  front  antenna 
(see  Figure  4.2),  where 

uD  =  2k  •  U  (3.18) 

and 

a  =  uD  At  =  y  ^  Ut  (3.19) 

where  B  is  the  baseline  separation  between  antennas,  V  is  the  aircraft  velocity,  Ur  is 
the  radial  c.  mponent  of  ocean  surface  motion  and  A  is  the  radar  wavelength.  Equa¬ 
tion  (3.19)  does  not  consider  speckle  contribution  from  the  random  ocean  surface 
and/or  from  the  measurement  and/or  processing  system  errors.  These  limitations 
are  reviewed  briefly  in  the  next  section  and  further  discussed  in  the  next  chapter. 

The  resulting  interferogram  is  a  complex  image  with  phase  which  is  the 
difference  of  the  two  original  SAR  image  phases  and  magnitude  which  is  the  product 
of  the  two  original  magnitudes.  Practically,  the  phase  and  amplitude  images  are 
processed  into  two  separate  images.  For  more  details  see  Chapter  IV.  B. 

2.  INSAR  Limitations 

The  limitations  of  INSAR  for  spatial  imaging  of  ocean  surface  velocity 
structure  may  be  categorized  into  two  types:  1)  system  dependent  2)  the  degree 
of  coherency  of  the  imaged  ocean  scene  within  the  inherent  time  measurement  of 
INSAR. 

a.  Platform  Stability 

In  order  to  compare  the  phases  of  a  rezel  (resolution  cell)  which  is 
observed  at  two  different  times,  the  viewing  aspects  must  be  nearly  the  same.  This 
might  not  be  fulfilled  due  to  lack  of  accurate  knowledge  of  relative  antenna  location, 
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which  in  turn  will  result  in  non-zero  velocities  and  acceleration  in  the  radial  direction 
(antenna  pattern  disimilarity).  For  Equation  (3.19)  to  extract  accurately  the  radial 
component  of  ocean  surface  motion,  it  is  assumed  that  the  two  antennas  are  collinear 
with  the  radar  motion.  Aircraft  pitch  and  yaw  will  add  phase  difference  between 
corresponding  pixels,  since  one  antenna  will  be  closer  in  range  to  the  observed 
resolution  cell.  The  radar  cannot  distinguish  between  a  slight  sideways  drift  of 
its  carrier  and  the  radial  motion  of  the  ocean  surface.  This  contamination  can  be 
compensated  using  readouts  of  aircraft  inertial  navigation  system  (INS),  or  to  a 
reasonable  degree,  by  inclusion  of  land  in  the  scene  to  give  absolute  reference,  or 
a  combination  of  the  two.  The  INS  provides  aircraft  pitch,  yaw,  pitch  and  yaw 
rates,  horizontal  and  vertical  velocity,  and  other  flight  parameters.  With  these 
data,  the  expected  error  phase  shifts  from  a  stationary  scene  can  be  computed  and 
compensated  accordingly.  A  yaw  error  of  0.01  degrees  will  yield  an  error  in  current 
estimates  of  4  cm  s-1  (Goldstein  and  Zebker,  1987). 

b.  Backscattering  Surface  Roughness  Coherence 

The  interferometric  technique  requires  comparison  of  two  samples 
of  a  moving  surface  element  (resolution  cell)  in  a  time  interval  A t  =  B/2V  required 
for  the  rear  antenna  to  move  (the  mean  aircraft  speed  is  220  m/s)  half  the  distance 
between  the  antennas  j3/2  =  9.9  m.  It  is  assumed  that  in  this  time  increment,  A i 
—  0.041  s,  the  surface  wave  field  as  a  whole  in  the  resolution  cell  (~  6  -r  8x  12  m) 
does  not  decay  or  grow,  i.e.,  the  scattering  patch  is  nearly  stationary.  A  number  of 
hydrodynamic  processes  affect  the  lifetime  of  Bragg  resonant  waves.  Among  these 
are  viscous  dissipation  and  their  response  (relaxation)  to  external  forcing,  such  as 
wind,  currents,  resonant  nonlinear  interactions,  straining  by  orbital  velocity  of  long 
waves,  and  turbulence.  The  viscous  half  life  time  decay  of  a  10  cm  wave  is  63  s.  But 
the  lifetime  of  wind  waves  are  controlled  by  relaxation  processes,  which  are  inversely 
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proportional  to  the  wind  growth  rate.  Therefore,  as  the  wind  speed  increases,  the 
relaxation  time  decreases.  For  10  cm  waves,  the  relaxation  time  is  about  two  seconds 
for  6  m/s  winds  (Valenzuela  and  Wright,  1976).  Hence,  for  moderate  sea  conditions 
the  ripples  should  remain  nearly  constant  within  the  time  interval  used  to  measure 
the  phase  difference.  However,  the  spread  of  orbital  velocities  in  the  radial  direction 
within  the  resolution  cell  will  also  decorrelate  the  backscatter  complex  signal  causing 
an  ambiguity  of  the  phase  as  the  ocean  wave  height  changes.  An  initial  attempt  to 
study  the  limitations  of  the  INSAR  technique  due  to  finite  scene  coherence  time  is 
made  in  the  next  chapter. 
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IV.  THEORETICAL  ANALYSIS  OF  INSAR  IN 
IMAGING  OCEAN  SURFACE  DYNAMICS 


A.  INTRODUCTION 

In  Chapter  III.  D,  the  known  physical  causes  for  the  ocean  surface  decorre¬ 
lation  were  described.  The  chapter  ended  with  a  recommendation  that  the  velocity 
spread  of  the  scattering  elements  due  to  local  orbital  velocity  of  the  subresolution 
scale  waves  should  be  further  analyzed.  In  this  section  for  the  first  time,  a  prelim¬ 
inary  theoretical  model  for  an  airborne  INSAR  imaging  the  ocean  surface  is  given. 
This  model  is  based  on  the  accepted  SAR  theory  of  imaging  the  moving,  nonrigid, 
distributed  ocean  surface  (e.g.,  Alpers  and  Rufenach,  1979;  Swift  and  Wilson,  1979; 
Valenzuela,  1980;  Raney,  1981;  Jain,  1981;  Rotheram,  1983;  Plant  and  Keller,  1983; 
Tucker,  1985;  Hasselmann  et  al.,  1985;  Monaldo  and  Lyzenga,  1986).  It  is  shown 
that  the  limited  ocean  surface  temporal  coherency  due  to  the  velocity  spread  of 
SAR  subresolution  waves  might  limit  the  capability  of  INSAR  in  providing  virtually 
instantaneous  high  resolution  snapshots  of  near  ocean  surface  velocity  patterns. 

In  the  analysis  given  below  the  SAR  two-scale  model  is  assumed  implicitly. 
This  model  was  introduced  first  by  Tucker  (1985)  who  modified  the  earlier  two-scale 
hydrodynamic  model  (Wright,  1968,  Valenzuela,  1978).  The  model  divides  the  ocean 
surface  waves  into  two  scales,  a  subresolution  short  scale  with  wave  numbers  greater 
than  the  SAR  azimuthal  cutoff  wavenumber  which  are  treated  statistically,  and  a 
long  scale  with  wave  numbers  shorter  than  this  wavenumber  which  are  treated  deter¬ 
ministically  (see  Chapter  III.  D).  It  is  found  that  the  motions  due  to  subresolution 
waves  produce  an  azimuthal  smearing  which  can  be  represented  by  a  Gaussian  low 
pass  filter  acting  on  the  azimuthal  component  of  wavenumber  in  the  image  (Tucker, 
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1985).  This  degradation  in  azimuth  resolution  (smearing)  is  inversely  proportional 
to  the  scene  coherence  time  and  determines  the  nominal  azimuth  resolution  for  SAR 
and  INSAR.  The  scene  coherence  time  for  L-band  radar  will  be  estimated  based  on 
the  Pierson-Moskowitz  ocean  spectrum  (350-600  msec)  and  the  backscatter  reflec¬ 
tivity  of  simultaneous  SAR  and  INSAR  observations  from  the  Marina  experiment 
on  September  8,  1989  (~  110  msec). 

A  major  accepted  simplification  of  SAR  imaging  theory,  both  for  mean  image 
and  speckle  statistics,  follows  from  the  hypothesis  that  the  complex  reflectivities  at 
different  points  on  the  sea  surface  are  uncorrelated.  The  spatial  separation  of  the 
different  scattering  elements  (facets)  should  be  large  in  comparison  with  the  Bragg 
wavelength,  but  still  small  in  comparison  with  the  SAR  resolution  scale.  The  surface 
reflectivity  may  therefore  be  regarded  as  spatially  white  (after  Hasselmann  et  al., 
1985).  This  assumption  will  be  used  in  the  dirivation  of  the  preliminary  model  of 
the  ocean  surface  cross-covariance  of  complex  reflectivity  for  INSAR. 

B.  INSAR  IMAGING  OF  THE  OCEAN 

The  complex  radar  return  signal  at  time  t  is  given  by  (Kasilingam  and 
Shemdin,  1988): 

S(t)=f  f  r(x,t)f(t  -2Ro/c)A(x  -Vt,y)dxdy  (4.1) 

J  — 00  J— oo 

with  x  =  x(a :,y),  see  Figure  4.1. 

where  r(x,  t )  is  the  time  variant  ocean  surface  Bragg  and  specular  reflectivity,  f(t  — 
2R0/c )  is  the  radar  pulse  function  and  A(x  -  Vt,y )  is  the  antenna  pattern.  We  shall 
ignore  the  compression  effect  in  the  pulse  function  since  the  time  variation  of  the 
returned  signal  within  the  pulse  duration  (few  microseconds)  is  negligible.  The  range 
dependence  of  the  pulse  function  introduces  a  quadratic  phase  term.  For  simplicity 
assume  that  the  pulse  rate  is  high  enough  to  fulfill  the  sampling  criteria  in  azimuth 
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and  a  continuous  signal  is  transmitted.  Therefore,  only  a  quadratic  phase  term  will 
be  introduced  by  the  pulse  function.  An  illustration  of  an  airborne  SAR  imaging 
the  ocean  surface  is  depicted  in  Figure  4.1  where  k  is  the  radar  wave  number  in 
the  viewing  direction,  V  is  the  platform  velocity  and  Rq  is  the  initial  aircraft  target 
slant  range  (see  Figure  4.1). 


Figure  4.1:  An  artistic  illustration  of  SAR  geometry  in  imaging  the  ocean. 

The  recorded  received  signals  are  processed  by  convolving  the  records  with  a 
matched  filter,  resulting  in  the  SAR  image.  For  a  focused  system  (see  Chapter  III) 
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with  no  focusing  adjustments  for  wave  pattern  motions,  the  matched  filter  impulse 
response  is  taken  to  be 

h(x  -  Vt)  =  exp  [tfc(x  -  Vt)2/Ro]  (4.2) 

The  complex  conventional  SAR  image  of  a  pixel  from  a  single  antenna  is 

» 

therefore  obtained  by  convolving  S(t )  given  by  Equation  (4.1)  with  the  matched 
filter  (Equation  (4.2)): 

/+00  f+oo  f+oo 

S(t  -  r)h(T)dr  —  /  /  r(x,t-r) 

•oo  J — OO  •/— oo 

•exp  |-ifc(x  -V(t-  r))2/Ro]  A(x  -V(t-  r)) 

•exp  [ifc(x  -  Vr)2/Ro]  dxdr  (4.3) 


where  the  first  exponential  term  on  the  right  hand  side  of  Equation  (4.3)  is  the 
quadratic  phase  history  for  a  stationary  target.  Similarly,  the  corresponding  complex 
pixel  image  for  the  second  antenna  taken  at  a  slightly  different  time  At  =  B/2V, 
which  is  the  time  required  for  the  rear  antenna  to  move  to  a  place  of  symmetry  with 
respect  to  the  forward  antenna,  is  given  by: 


/+oo 

S{t  +  B/2V  -  rXrOdr' 
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•  exp  [ik{x  -  V{t  +  B/2V  -  t'))2/R ]  A{x  -  V{t  +  B/2V  -  t')) 
•exp  \ik(x  -  Vt')2 / dx1  dr'  (4.4) 


Equations  (4.3),  and  (4.4)  represent  complex  SAR  images  of  the  ocean  surface 
reflectivity  of  the  same  pixel.  The  goal  is  to  extract  reliably  the  mean  Doppler  shift 
caused  by  the  first  order  ocean  surface  motion.  This  phase  information  is  included 
entirely  in  the  cross-covariance  of  the  complex  surface  radar  reflectivity  (r^)  (see 
Equation  (4.6)  below).  From  a  statistical  analysis  point  of  view,  the  averages  of 
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rir2  (denoted  by  brackets)  are  defined  for  an  infinite  ensemble  of  different  surfaces. 
Practically,  the  averages  are  estimated  over  a  finite  surface  region  for  a  finite  time. 
By  taking  the  cross-covariance  of  the  corresponding  complex  images  from  the  two 
antennas,  the  phase  difference,  which  is  proportional  to  the  virtual  surface  velocities, 
is  obtained.  The  interferometric  operation  is  therefore  given  by: 

Q(x)  =  (n(t  =  x/V)  i;(t  +  B/2V))  (4.5) 

Substituting  Equations  (4.3),  (4.4)  into  (4.5)  yields: 

/+00  f+OO  f+OO  f+OO 

III  (ri{x,t-r)r;(x',t  + B/2V -t')) 

-oo  J—oo  00  J  —  oo 

•exp  [— ik(x  —  V(t  —  T))2/i2o]  Mx  —  ~  T ))  exP  \ik(x  —  Vr)2 /ilo] 

•exp  [ik{x'  -  V{i  +  B/2V  -  t'))2 / R^Aix'  -  V(t  +  B/2V  -  t')) 

•exp  ^ik(x' —VT,)2/Ro)dxdxt drdr'  (4.6) 

The  first  term  on  the  right  hand  side  of  Equation  (4.6)  is  the  covariance  of 
complex  reflectivity  of  the  ocean  surface,  which  is  defined  by  the  second  moment 
statistics  of  the  ocean  surface  backscattered  signal.  As  noted  in  the  introduction 
to  this  section,  the  averaged  SAR  mean  image  is  directly  related  to  these  statistics, 
namely,  the  covariance  of  complex  reflectivity  of  the  ocean  surface.  Figure  4.2 
illustrates  the  geometry  of  INSAR. 

The  extraction  of  the  phase  information  from  the  backscattered  moving  sur¬ 
face  by  using  Equation  (4.5)  is  the  essence  of  INSAR  in  imaging  the  ocean.  This  is 
fundamentally  different  from  the  SAR  image  intensity,  given  by: 

<w2>  =  <*'i*;>  m 

Equation  (4.7)  clearly  illustrates  the  loss  of  the  phase  information  in  the  SAR  image, 
which  is  retained  in  Equations  (4.5),  (4.6)  for  INSAR. 
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Figure  4.2:  Geometry  of  INSAR:  a)  top  view,  b)  perspective  view  shows 
the  change  in  distance  of  the  resolution  element  due  to  wave  motion  in 
the  time  it  takes  for  the  rear  antenna  to  more  to  a  symmetric  position 
with  respect  to  the  front  antenna. 

The  SAR  System  Transfer  Function  (STF)  is  not  considered  at  the  present 
stage.  It  is  assumed  that  the  SAR  is  a  perfect  instrument  that  maps  the  covariance 
of  complex  reflectivity  into  a  two-dimensional  complex  image.  The  interest  here  is  to 
investigate  the  feasibility  of  INSAR  in  imaging  the  moving  ocean  surface.  Therefore, 
the  interferometric  operation  described  by  Equation  (4.5)  can  be  replaced  by  the 
cross-covariance  of  complex  reflectivity  of  an  imaged  ocean  scene  sensed  by  the  two 
INSAR  antennas  with  a  time  difference  A t  =  Bj 2V. 

C.  DERIVATION  OF  THE  CROSS  COVARIANCE  OF 
THE  COMPLEX  RADAR  REFLECTIVITY  FOR 
INSAR 

The  complex  radar  reflectivity  for  an  individual  scattering  element  moving 
with  velocity  U(x,t)  is  given  by: 

r(x,  t )  =  r0(x,  t)  exp  [2ik  •  U(x,  t )  (4.8) 
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where  ro(x,  t)  is  the  instantaneous  complex  reflectivity  of  an  ocean  surface  scattering 
element  within  the  SAR  resolution  cell.  As  noted  earlier,  these  scattering  elements 
(facets)  are  assumed  to  be  large  relatively  to  the  resonant  backscattering  Bragg 
waves  (5-10  times),  but  much  smaller  than  the  SAR  resolution  cell  (rezel). 

The  cross-covariance  of  complex  reflectivity  of  a  time  dependent  backscat¬ 
tering  surface  f::om  corresponding  rezels  as  imaged  by  the  two  antennas  is  given 
by: 

(ri(x, t)  r^x',*'))  =  (r0i(x,t)exp  [2tk-U(x,<)  t]  r^x'.t') 

•  exp  [-2ik  •  U(x',  t1)  f'j  )S(x  -  x')  cr0(x,  t )  (4.9) 

where  t'  —  t  +  At  =  i  +  B/2V  and  6(x  —  x')  represents  the  white  spatial  distribution 
of  the  scattering  elements  within  the  SAR  resolution  cell  (Hasselmann  et  al.,  1985). 
<r0(x,/)  is  the  variance  of  r  and  represents  the  modulation  of  the  covariance  of 
short  wave  reflectivity  by  the  slow  time  varying  long  wave  tilt  and  hydrodynamic 
modulation.  oo(x,  t)  represents  the  radar  cross-section  as  imaged  by  conventional 
SAR  and  real  aperture  radars. 

Ir  Equation  (4.9),  the  l  me  lag  A t  represents  only  the  constant  time  interval 
required  for  INSAR  measurement.  The  additional  time  difference  of  r  —  r'  (see 
Equation  (4.6)),  which  is  of  the  order  of  the  SAR  integration  time,  is  not  considered 
here  explicitly  (except  when  determining  the  effective  azimuth  resolution).  This 
time  difference  is  not  relevant  for  INSAR  operation  once  the  two  corresponding 
complex  SAR  images  are  formed.  For  more  detail  on  the  covariance  function  for 
SAR,  see  for  instance  Hasselmann  et  al.  (1985).  The  discussion  here  (Equation 
(4.9))  is  based  on  analyzing  the  second  order  statistics  of  the  fluctuating  reflectivity 
represented  by  the  ensemble  term  in  Equation  (4.9)  on  the  time  scale  relevant  for 
INSAR  {M  =  B/2V). 
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The  SAR  two-scale  model  is  invoked  resulting  in  the  ocean  surface  velocities 
decomposed  as 

U(x,  t)  =  Uo(x,  t)  +  U'(x,  t)  (4.10) 

where  U0(x,  t)  is  the  long  scale  velocity  and  U'(x,  t)  is  the  subresolution  short  scale 
velocity.  Uo(x,f)  can  be  extracted  by  INSAR  if  the  temporal  variability  of  the 
backscattering  surface,  defined  by  the  radar  scene  coherence  time,  is  small  relative 
to  the  time  difference  At  =  B/2V  required  to  obtain  the  phase  image  of  INSAR 
which  is  directly  proportional  to  the  mean  radial  surface  velocity. 

The  limits  of  the  temporal  decorrelation  of  the  ocean  surface  on  INSAR  in 
obtaining  the  surface  velocity  Uo(x,  f),  are  now  addressed.  The  variation  of  the 
ocean  surface  velocity  U(x,f)  in  the  time  At  =  B/2V  can  be  described  by  using 
the  Taylor  expansion 

U(x,i  +  B/2V)  =  U(x,f)  +  B/2V  (4.11) 

Substituting  Equations  (4.10)  and  (4.11)  into  Equation  (4.9)  and  factoring  out  the 
deterministic  long  scale  contributions  yields 


(ri  (x,f)  rj(x,i  +  B/2V)) 


exp[-2zk  •  U0(x,f)  B/ 2Vr] 

•exp[— 2ik-^  (x,t)  ~  ( +  <)]^o(x,f) 


roi(x,t)  rj2(x,t  +  B/2V) 


(4.12) 


•exp 


_2i— k 
2V 


U'(x,*)  + 


dU'(x,f) 

dt 


t  + 


B 

2V 


Equation  (4.12)  represents  the  fundamental  basis  of  the  present  study  in  evaluating 
the  feasibility  of  INSAR  in  imaging  the  ocean  surface.  It  combines  the  relevant 
physics  of  the  backscattering  surface  as  a  function  of  INSAR  physical  configuration. 
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The  second  exponential  term  on  the  right  hand  side  of  Equation  (4.12)  can 
be  neglected  if 

<4-13> 

where  t  is  of  the  order  of  the  SAR  integration  time.  Since  t  »  B/2V  and  the 
absolute  value  of  the  acceleration  of  the  dominant  wave  is  giver  by  0Uo /dt  —  wllo, 
inequality  (4.13)  can  be  rewritten  as: 

ut  <  1  (4.14) 

or 

T  »  2 irt  (4.15) 

where  T  is  the  period  of  the  imaged  dominant  wave.  Equation  (4.15)  specifies  the 
limiting  condition  for  the  integration  time  when  imaging  dominant  waves  by  SAR 
and  INSAR.  The  effective  processing  integration  time  for  mean  incidence  angles 
(30°)  is  almost  1.7  sec.  In  practice,  the  scene  coherence  for  a  given  radar  wavelength 
dictates  the  effective  integration  time,  which  is  substantially  less  than  1  sec  (e.g., 
Tucker,  1985).  This  crucial  issue  is  discussed  in  the  next  two  sections. 

The  first  term  on  the  right  hand  side  of  Equation  (4.12)  is  the  required  phase 
shift  to  be  determined  by  INSAR.  It  represents  the  first  order  ocean  surface  radial 
velocities.  U0  is  a  vector  sum  of 

Uo  =  Ucp  +  Uob  +  Uch  +  Var  (4.16) 

The  individual  contributions  are  identified  as: 

•  Ucp  phase  velocity  of  the  backscattering  resonant  Bragg  waves; 

#  Uob  mean  orbital  velocities  of  waves  substantially  longer  than  the  character¬ 
istic  dimension  of  the  resolution  cell. 
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•  Uc/i  material  velocity  of  any  near  surface  currents  underlying  the  backscat- 
tering  Bragg  waves. 

•  Uytfl  aircraft  attitude  error  velocity  due  to  yaw  and  pitch  bias  of  the  INSAR 
antennas  (see  Chapter  III.  B). 

Practically,  the  ensemble  average  factor  in  Equation  (4.12)  expresses  the  ad¬ 
ditional  fading  in  the  cross  covariance  of  complex  reflectivity  due  to  the  slight  time 
difference  At  when  the  same  scene  is  imaged  separately  with  a  small  time  delay 
by  the  two  INSAR  antennas.  This  phenomena  is  defined  here  as  c(x,  t,  At).  These 
modifications  in  Equations  (4.13),  (4.15)  substituted  into  Equation  (4.12)  yield: 

(ri(x,t)  r2(x,<  +  B/2V))  =  exp[-2tk-U0(x,<)  B/2V]c(x,i,  B /2V)ao(x,t)  (4.17) 

Both  are  shown  to  be  quantified  as  a  function  of  the  finite  scene  coherence 
time.  c(x,t,  At)  can  be  rewritten  as  (see  Equation  (4.12)  (4.17)): 

c(x.t,B/2V)  =  exp  ~2i~k  •  fu'(x,t)  +  (1  +  ^))  c'(x> *< B/2V) 

(4.18) 

with 

c'(x,  t,  B/2V)  =  (r0I(x,  t)  r*2(x,  t  +  B/2V))  (4.19) 

where  d(x,t,B/2V)  is  the  temporal  cross- correlation  function  value  for  the  INSAR 
time  difference  At  =  B/2V.  The  phase  term  in  Equation  (4.17)  can  be  rewritten  in 
terms  of  radian  Doppler  frequency 

vls  =  2k  •  U'(x,t)  (4.20) 

where  u3S  is  the  perturbation  in  radian  Doppler  shift  associated  with  the  velocity 
spread  of  modulating  the  subresolution  scale  wave  scattering  elements.  Substituting 
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Equation  (4.20)  into  Equation  (4.18)  yields: 


c(x,t,B/2V)  =  exp 


.  B 
l2V 


(mo +*#&*) 


c'(x,t,B/2V)  (4.21) 


where  dut,fdt  ( B/2V )2,  which  is  O(10-2),  is  neglected. 

In  contrast  with  SAR  image  intensity  which  ignores  phase  information  (Equa¬ 
tion  (4.7)),  INSAR  imagery  represents  the  virtual  phase  difference  between  the  two 
corresponding  SAR  images.  Therefore,  the  present  model  (Equations  (4.14-4.21)) 
takes  into  account  additional  complex  terms  which  were  not  considered  explicitly  in 
SAR  models.  We  shall  further  concentrate  in  evaluating  the  complex  disturbance 
to  INSAR  in  sensing  ocean  surface  velocities  introduced  by  c(x,  t,  At).  This  dis¬ 
turbance  includes  cross  section  attenuation  and  additive  phase  noise,  both  caused 
by  the  subresolution  intermediate  waves  modulating  randomly  the  backscattering 
elements  (facets).  The  finite  temporal  coherence  of  the  ocean  scene  is  of  special  im¬ 
portance  since  it  defines  the  maximum  coherent  radar  integration  time  for  a  given 
radar  wavelength,  and  therefore  the  maximum  azimuth  resolution  achievable.  In  the 
p.  esent  study,  the  effect  of  the  finite  scene  coherence  in  the  time  interval  (B/2V) 
required  to  obtain  the  interferometric  information  based  on  the  corresponding  inte¬ 
grated  SAR  images  from  the  two  antennas  is  discussed. 

When  considering  the  radar  scene  coherence  time,  two  time  scales  should  be 
looked  at: 


•  The  decorrelation  time  of  an  individual  facet  within  the  SAR  resolution  cell. 
This  fading  in  the  mean  complex  reflectivity  can  be  described  by  a  temporal 
decorrelation  function. 

•  Intrinsic  coherence  time  in  which  the  Bragg  waves  of  a  single  scattering  el¬ 
ement  (facet)  preserve  their  structure.  It  is  assumed  that  this  time  scale  is 
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long  compared  with  the  coherent  integration  time  (Hass.,  1985,  Vachon  et  ah, 
1988). 


The  temporal  correlation  function  is  the  inverse  Fourier  transform  of  the 
reflectivity  variance  Doppler  spectrum.  Based  on  theory  (e.g.,  Tucker,  1985)  and 
the  TOWARD  experiment  (1984),  it  can  be  assumed  that  this  Doppler  broadening 
of  short  wave  reflectivity  has  a  Gaussian  shape.  After  removing  the  mean  Doppler 
shift  due  to  long  wave  and  current  velocities,  the  velocity  power  spectrum  of  the 
short  wave  reflectivity  from  a  given  pixel  is 

■P(w)  =  J  (r0i(t)  r*02(t  -f  r))  exp [-iwrjdr  =  P0exp  [-u£/2a*] ,  (4.22) 

where  aw  is  half  Doppler  bandwidth  of  P(w). 

The  inverse  Fourier  transform  of  Equation  (4.22)  yields  the  temporal  cross¬ 
correlation  function  of  short  wave  reflectivity  at  a  given  point. 

1  r+°°  r  i 

c'(t)  =  —  /  P( u>)exp  uwrldu;  (4.23) 

Substituting  Equation  (4.22)  into  Equation  (4.23)  yields 

C#(T)  =  ^  /  exP  [-w2/2c£]  exp  t’wrjcfu;,  (4.24) 

where  u)  =  uSi. 

Performing  the  integration  in  Equation  (4.24)  renders 

c'(r)  =  JL/Etll  exp  \-vyn]  =  ^  exp  [-«7^2/2]  (4.25) 

The  standard  deviation  au  of  the  Doppler  broadening  of  short  wave  reflec¬ 
tivity  represents  the  inverse  of  the  correlation  time  of  the  subresolution  scattering 
elements  (facets),  namely,  the  finite  scene  coherence  time. 


(4.26) 
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Substituting  Equations  (4.25)  and  (4.26)  into  Equation  (4.21)  yields 
c{x,t,r)  =  Ph  exp  exp  [tr  (u„  +  , 


c(x,  t,  r)  =  PI)  exp  p exP  |^T  +  ~^‘tj  j ’  (4-27) 

where  Pq  =  <juiPo(2w)~1^2  and  r  =  B/2V. 

The  phase  term  in  Equation  (4.27)  expresses  the  Doppler  spread  of  the  scat¬ 
tering  element  (facet)  bandwidth,  which  is  assumed  to  be  the  same  for  all  facets. 


Thus 


+ 


du..  t 


0(  Aw)  = 


Substituting  Equation  (4.28)  into  Equation  (4.27)  yields 


,  ,  v  D/  T2]  \i2nr 

c(x,<,t)  «  Pq  ex p  —  —  exp  - 

2tcJ  rc 


(4.28) 


(4.29) 


Equation  (4.29)  represents  explicitly  the  previously  mentioned  complex  dis¬ 
turbance  to  INSAR  remotely  sensing  the  ocean  surface.  Substituting  Equation 
(4.29)  into  Equation  (4.17)  yields  the  explicit  expression  of  the  cross-covariance  of 
complex  reflectivity  for  SAR  interferometry. 

(r^x,  t)  rj(x,<,  At))  =s  <r(x,  t)  exp  — 

rc 

•exp  -—p  exp  [— 2ik  •  U0(x,  f )  At]  (4.30) 


"eXP[  2r2j  exPl-ZiK’  uolx>tJ  Atj 

where  At  =  B/2V.  Equation  (4.30)  illustrates  the  sensitivity  in  measuring  the  line 
of  sight  velocity  component  U0  of  an  unsteady  rough  ocean  surface.  It  is  shown 
that  the  key  limiting  parameter  is  the  ratio  between  the  inherent  interferometric 
time  delay  (At  =  B/2V )  and  the  scene  coherence  time  rc  which  is  analyzed  next. 

The  first  phase  term  in  Equation  (4.30)  represents  the  additive  phase  noise  to 
the  mean  surface  velocity  Uo,  which  is  attempted  to  be  measured  by  INSAR.  This 
noise  term  is  required  to  be  small  relative  to  the  mean  phase  term.  This  implies 


|2k  •  U0At|  » 


I  27T  At  I 


(4.31) 


which  results  in 


Tc>2|U0|' 

For  the  NASA/JPL  DC-8  airborne  configuration,  the  range  of  measurable  velocities 
Uo  is  ±1.25  m/s  and  for  the  L-band  radar  wavelength  of  A  =  24  cm,  Equation 
(4.32)  indicates  that  the  scene  coherence  time  rc  should  be  at  least  0.1  sec  for  the 
additive  phase  noise  term  to  be  neglected.  The  known  observed  estimates  of  the 
scene  coherence  time  for  L-band  radar  are  equal  and  smaller  than  the  order  of  102 
msec  (Lyzenga  and  Shuchman,  1983,  Kasilingam  and  Shemdin,  1988). 

Next  to  be  examined  is  the  attenuation  of  the  radar  cross-section  caused  by 
the  slightly  different  time  in  imaging  a  given  scene  by  the  two  physically  separated 
INSAR  antennas.  This  attenuation  is  represented  by  the  second  exponential  in 
Equation  (4.30).  The  radar  cross-section  for  INSAR  can  be  expressed  as: 


</(x,  f ,  At)  =  <ro(x,  t)  exp(-A/2/2r2)  (4.33) 

In  order  to  extract  the  phase  information  representing  the  surface  velocity 
(last  term  on  RHS  of  Equation  (4.30)),  one  needs  a  sufficient  signal  <r'(x,t,At) 
above  the  noise  level.  A  3  dB  attenuation  of  the  cross-covariance  ‘signal’  is  allowed, 
ensuring  accurate  measurements  of  the  backscattering  surface  velocity.  The  required 
scene  coherence  time  rc  in  Equation  (4.30)  is  taken  to  be 

Tc  ~  2V  (2  In  2)1/2  (4‘34^ 

Once  tc  is  known,  the  upper  limit  on  At  =  B/2V  is  defined,  which  in  turn  prescribes 
the  necessary  INSAR  set  up.  Conversely,  for  a  given  INSAR  set  up  the  minimum 
allowable  scene  coherence  time  is  determined.  For  the  DC-8  INSAR  configuration 
( B  =  19.8  m,  platform  velocity  V  «  200  m/s),  using  criterion  (4.34)  yields  the 
minimum  required  scene  coherence  time  rc  =  20  msec.  Combining  the  phase  and 
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amplitude  constraints  on  the  radar  scene  coherence  time,  the  temporal  correlation 
of  the  short  wave  reflectivity  should  be  substantially  greater  than  100  msec  for  the 
given  INSAR  system. 

Finally,  the  model  can  estimate  the  range  of  variability  of  the  complex  re¬ 
flectivity  within  the  inherent  INSAR  measurement  time,  which  will  still  allow  the 
extraction  of  ocean  surface  resolvable  motion  with  minimum  contamination. 

D.  ESTIMATION  OF  SCENE  COHERENCE  TIME 

The  next  step  in  the  study  is  to  estimate  the  actual  scene  coherence  time  for 
an  L-band  imaging  radar  using  two  different  approaches: 

•  Theoretically,  assuming  the  Pierson  Moskowitz  ocean  spectrum. 

•  Experimentally,  using  for  the  first  time  simultaneously  SAR  and  INSAR  imag¬ 
ing  records. 

1.  Estimation  of  Scene  Coherence  Time  Based  on 
Pierson-Moskowitz  Spectrum 

The  finite  scene  coherence  time  may  be  one  of  the  most  important  factors 
in  constraining  SAR’s  ability  to  image  azimuthally  traveling  waves  (see  Section  IV. 
B).  In  this  section,  the  analysis  will  follow  recent  attempts  to  quantify  this  parameter 
(see  e.g.,  Tucker,  1985b,  Alpers  and  Bruening,  1986,  Vachon  et  al.,  1988)  using  the 
Pierson-Moskowitz  ocean  spectrum  for  fully  arisen  sea.  Kasilingam  and  Shemdin 
(1988)  point  out  the  uncertainty  in  the  validity  of  this  approach,  which  might  lead  to 
differences  by  order  of  magnitudes  in  estimating  the  scene  coherence  time  rc.  Due  to 
the  complexity  of  the  temporal  ocean  surface  decorrelation  rc,  the  theory  associated 
with  the  ocean  scattering  motion  effects  is  still  under  extensive  research.  Owing  to 
the  lack  of  experimental  evidence  and  more  profound  theoretical  understanding  of 
the  small  scale  roughness  of  the  sea  surface,  this  approach  still  remains  acceptable. 
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The  classical  approach  is  to  relate  the  azimuth  image  smear 
{6{x)2)'t*  (which  is  inversely  proportional  to  the  scene  coherence  time),  to  the  rms 
facet  surface  velocities  (U'2)1/2  caused  by  the  random  material  velocities  of  subres¬ 
olution  scale  waves  (see  Equation  4.20). 

(6(x)2yt2  =  (u  yy^R/v  (4.35) 

where  (I/'2)1/2  is  the  radial  component  of  (U ,2)1/2. 

The  shape  of  the  Doppler  broadening  spectrum  associated  with  the  ran¬ 
dom  facet  velocities  U'  is  assumed  to  be  Gaussian  (see  Equation  (4.22)),  resulting 
in  a  SAR  azimuthal  imaging  transfer  function  with  a  low-pass  filter  character.  The 
bandwidth  of  this  filter  is  inversely  proportional  to  the  decorrelation  time  scale, 
which  is  interpreted  as  the  scene  coherence  time.  In  the  wave  number  domain,  the 
width  of  this  filter  (transfer  function)  is  inversely  proportional  to  the  azimuthal  im¬ 
age  smear.  It  is  important  to  note  that  the  loss  of  resolution  due  to  finite  scene 
coherence  time  is  not  a  defocusing  effect,  and  is  therefore  not  recoverable  by  focus 
adjustment  in  the  processing  or  integration  time  setting  (see  Chapter  III). 

Based  on  the  relation  between  the  scene  coherence  time  and  rms  surface 
velocity  (U'2)1/2  spread,  the  subresolution  velocity  spread  is  first  estimated.  (U'2)1/2 
may  be  computed  by  integrating  the  orbital  velocity  power  spectrum  between  the 
appropriate  wave  frequencies 

(U;2)  =  fw25»  /2(7)  du,  (4.36) 

•'W 1 

where  5„(w)  is  the  omnidirectional  surface  elevation  frequency  power  spectrum  and 
/(7)  is  the  geometric  function 

f(n)  =  —  f  (sin2  7  sin2<6  +  cos2  (4.37) 

7T  J-ir/ 2 
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where  <j>  is  the  angle  between  the  wave  train  and  the  aircraft  attitude  and  7  is  the 
incidence  viewing  angle  (Figure  4.1).  The  geometric  function  7(7)  represents  the 
projection  of  the  moving  scene  in  the  radial  viewing  direction.  7(7)  is  always  smaller 
than  1. 

A  rather  crude  assumption  is  next  made  about  the  ocean  wave  spectrum 
in  calculating  (U'2)1/2,  which  is  in  fact  sensitive  to  local  sea  conditions.  In  the 
absence  of  relevant  experimental  data,  it  is  assumed  that  the  ocean  wave  spectrum 
for  a  fully  arisen  sea  can  be  described  by  a  Pierson- Moskowitz  spectrum  (e.g.,  Alpers 
and  Bruening,  1986,  Tucker,  1985b).  Substituting  the  Pierson- Moskowitz  spectrum 
into  Equation  (4.36)  yields: 

(U'2)=jf  w2ag2u~sexp  [— /%/Uw)4]  f2{^)duj  (4.38) 

where  a  =  8.1  x  10-3  is  the  Phillips  constant,  =  0.74,  and  Uw  is  the  wind  speed 
at  a  height  of  19.5  m. 

After  integration, 

/  \  l!2 

{V?)  =  j  (I)  US  erf  [„  (U.)i  (4.39) 

with  n  =  (PY^2g2  =  82.78  m2  sec"2. 

The  lower  limit  of  integration  is  associated  with  the  shortest  resolvable 
ocean  wavelength  0  (10  m)  and  is  determined  by  the  Nyquist  criterion.  The  upper 
limit  is  determined  by  the  fundamental  subresolution  scattering  elements  (facet) 
with  size  of  the  order  of  Bragg  wavelength  (~  0.2  m)  which  can  be  replaced  for 
computational  convenience  by  w2  — »  oo  with  a  small  error. 

The  lower  limit  is  controversial  since  the  actual  SAR.  (and  INSAR)  az¬ 
imuth  resolution  is  not  known  a  priori,  being  a  function  of  the  scene  coherence 
time.  A  few  choices  based  on  different  physical  considerations  will  be  discussed 
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below.  The  lower  integration  limit  u>i  in  Equation  (4.39)  is  replaced  by  a  cutoff 
resolvable  wavelength  L\  using  the  dispersion  relation  for  water  waves, 

u>i  =  (2ng/L1)1/2  (4.40) 


Since  the  shortest  resolvable  wavelength  is  0  (10  m),  the  surface  tension  effects  are 
neglected  and  only  gravity  waves  are  considered.  By  inserting  the  integration  limits 
into  Equation  (4.39),  it  is  found 


(4.41) 


with  n'  =  n/g2ir  =  1.35  m  sec-2.  For  strong  winds  /i'Xi/U2  <  1,  the  radial  velocity 
spread  from  Equation  (4.41)  simplifies  to: 


(u;2)1'2 


a 

4 


m 


(4.42) 


which  is  independent  of  the  sea-state.  For  weak  wind  speeds  y!Lx/ U2  >  1  and 
Equation  (4.41)  simplifies  to: 


a 

4 


Uw  /( 7) 


(4.43) 


Surprisingly,  the  rms  velocity  spread  for  weak  winds,  Equation  (4.43),  is  indepen¬ 
dent  of  the  cutoff  wavelength  Lx,  which  is  determined  in  this  case  mainly  by  the 
SAR  system.  This  wavelength  represents  the  scale  limit  of  waves  which  are  treated 
deterministically  ( L  >  L\)  and  waves  which  are  treated  statistically  (L  <  L\).  This 
statistical  region  is  of  interest  since  the  unresolved  wave  (L  <  L\)  modulate  ran¬ 
domly  the  fundamental  scattering  elements.  The  temporal  variability  in  the  surface 
reflectivity  is  determined  by  the  scene  coherence  time. 

Two  related  choices  for  the  cutoff  wavelength  L\  have  been  considered. 
The  first  choice  is  based  upon  the  processor  coherent  integration  time  for  a  wave 


50 


with  a  period  equal  to  T.  This  wave  and  shorter  waves  will  propagate  at  least 
one  cycle  over  the  time  of  image  formation  contributing  only  to  the  rms  velocity. 
Using  again  the  deep  water  dispersion  relation  for  gravity  waves,  this  cutoff  scale  is 
(Vachon  et  al.,  1988): 

Li  =  £rl  (4.44) 

The  second  choice  used  by  some  investigators  (e.g.,  Alpers  and  Bruening,  1986, 
Tucker,  1985a)  is  L\  =  2 p,  twice  the  matched  filter  processor  resolution.  For  a  mov¬ 
ing  ocean  surface,  the  azimuth  resolution  is  degraded  and  is  practically  unknown. 
It  could  be  much  larger  than  2 /?,  which  is  the  nominal  azimuth  resolution  for  an 
unmoving  scene.  Noting  that  the  mean  Doppler  broadening  due  to  the  velocity  of 


the  intermediate  unresolvable  waves  is  given  by: 


fdc  = 


2<u;2)i/2 


(4.45) 


leads  to  the  finite  radar  scene  coherence  time 

Tc  =  7-.  (4.46) 

Jdc 

For  strong  wind  (Uw  >  10  m/s)  substituting  Equation  (4.44)  into  Equa¬ 
tion  (4.39)  and  the  result  into  (4.45)  and  (4.46),  yields  the  scene  coherence  time. 


t.  =  5.04 


77(7) 


and  for  Lx  =  2 p 


tc  =  4.45 


(4.47) 


(4.48) 


C  (pY'2f(  7) 

For  weak  winds  (Uw  <  5  m/s)  substituting  Equation  (4.43)  into  Equa¬ 
tion  (4.45),  (4.46)  gives: 


rc  =  7.74 


Uw  /(t) 


(4.49) 


TABLE  4.1:  Estimates  for  Scene  Coherence  Time  as  a  Function  of  Wind 
Speed  and  Cutoff  Wavelength  for  NASA/JPL  DC-8  L-Band  Radar  with 
/( 7)  =  1 


Equation 

lx 

Uw  (m/sec) 

rc  (msec) 

(4.47)1 

J-T2 

>  10 

600 

(4.48)2 

2  P 

>  10 

610 

(4.51)3 

350 

1  T  =  2  sec  2  p  =  3.03  m  3  H,  =  0.6  m. 


Since  for  fully  developed  wind  seas,  the  wind  speed  Uw  at  the  height  of  19.5  m  is 
related  to  the  significant  waveheight  (inks  units)  by 

=  47.4ff„  (4.50) 

Equation  (4.49)  can  be  rewritten  as 

Tt = hl2w. wr  (4'51) 

The  estimated  values  of  the  scene  coherence  time  in  Table  4.1  using 
Equations  (4.47),  (4.48),  (4.51)  and  /(j)  =  1  are  rather  long  compared  to  the 
sparse  estimates  based  on  observations  in  the  literature  (30-100  msec  for  L-band 
radar).  The  decrease  in  /( 7)  with  increasing  incidence  angle  will  increase  even  more 
the  scene  coherence  time  rc  as  computed  by  the  equations  used  in  Table  4.1.  The 
dependence  of  tc  on  the  spreading  function  7(7)  is  illustrated  in  Figure  4.3.  This 
dependence  was  not  considered  explicitly  in  other  works  and  its  importance  rises 
from  analysis  of  observed  data  (see  Chapter  VII). 

Variability  of  several  hundred  msec  in  the  estimates  of  coherence  time 
scales  is  noted  in  Table  4.1.  The  first  two  estimates  are  computed  using  the  DC-8 


52 


Figure  4.3:  The  effect  of  spreading  function  f(~f)  on  the  scene  coherence 
time  rc,  showing  the  increase  of  the  scene  coherence  time  with  incidence 
angle. 

SAR  parameters  and  are  independent  of  sea  state,  a  situation  which  seems  to  be 
unlikely.  The  third  estimate  seems  to  be  more  sensitive  to  the  real  ocean  situation, 
namely,  dependence  of  scene  coherence  time  on  sea  state  and  hence  wind  speed.  It 
will  be  shown  through  the  experimental  approach  adopted  in  the  next  section  that 
the  shortest  coherence  time  estimate  of  the  present  method  is  rather  longer  than 
that  which  is  observed.  This  might  indicate  the  inadequacy  of  the  present  approach 
to  estimate  the  scene  coherence  time,  mainly  because  of  a  lack  of  comprehensive 
understanding  of  short  wave  radar  reflectivity  (Hasselmann  et  al.,  1985).  To  have  a 
more  reliable  estimate  of  this  parameter,  it  is  crucial  to  obtain  more  direct  measure- 
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merits  of  the  coherence  time  which  determines  the  effective  azimuthal  resolution  of 
the  SAR  and  hence  INSAR  imagery. 

2.  Estimation  of  Radar  Scene  Coherence  Time  Based  on  SAR  and 
INSAR  Simultaneous  Images 

The  temporal  coherence  of  the  ocean  scene  can  be  investigated  by  using 
simultaneous  SAR  and  INSAR  records.  As  described  in  section  B  of  this  chapter, 
the  image  intensity  of  SAR  is  a  map  of  the  covariance  of  the  complex  reflectivity. 
The  mean  image  intensity  obtained  by  a  single  antenna  is  computed  by  ensemble 
averaging  neighboring  pixels,  and  is  directly  related  to  the  covariance  of  complex 
reflectivity: 

{rx(x,<K(x,t  +  r))  =  <7i(x,t)exp  [-r2/2rc2]  (4.52) 

Similarly,  for  the  second  antenna  separated  physically  by  a  distance  B  from  the  first 
antenna,  the  covariance  of  complex  reflectivity  is  taken  to  be: 


(r2(x,  t  -f  At)r!J(x,  t  +  At  +  r))  =  <r2(x,  t ,  At)  exp 


(4.53) 


where  At  =  B/2V  and  cr  is  the  radar  cross  section  which  represent  the  covariance 
of  short  wave  reflectivity  modulated  mainly  by  long  wave  tilt.  Since  the  temporal 
variations  of  these  waves  are  slow  compared  to  the  interferometric  time  interval  At, 
it  is  reasonable  to  assume  that  the  radar  cross-section  does  not  vary  during  the  time 
At.  Therefore 

oi(x,f)  ~  cr2(x,t  +  At)  ~  a(x,t)  (4.54) 

When  corresponding  pixels  from  the  two  images  are  combined  to  extract 
the  interference  phase,  the  temporal  cross-correlation  between  the  corresponding 
pixel  scattering  elements  introduces  an  additional  cross-section  modulation  due  to 
the  random  velocity  spread  of  subresolution  intermediate  waves  modulating  the 
short  wave  reflectivity.  This  relatively  fast  time  decorrelation  of  the  scene  in  a  given 
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time  interval  At  =  B/2V  attenuates  the  INSAR  signal  (see  Equation  {4.33)).  The 
magnitude  of  the  cross-product  of  the  complex  images  is  directly  related  to  the 
cross-covariance  of  complex  reflectivity,  and  is  given  by: 

(ri(x,i)r;(x,t  +  At  +  r))  =  <rll2(x,t  +  At)  exp  [-r2/2rc2]  (4.55) 

where 

o-it2(x,  t,  At)  =  cr(x,  i)  exp  [- At2/2r2]  (4.56) 

Substituting  Equation  (4.56)  into  Equation  (4.55)  yields: 

(ri(x,t)rl(x,t  +  At  +  r))  =  a(x,i)exp  [-At2/2r2]  exp  [-r2/2r2]  (4.57) 

Neglecting  for  a  moment  the  additive  noise  effects,  the  effective  radar 
cross-section  obtained  by  INSAR  (Equation  4.56)  normalized  to  the  product  of  the 
square  root  of  both  corresponding  SAR  radar  cross-section  by  invoking  Equation 
(4.54)  yields 

K  -  (ri(x,t)r^(x,t  +  At)) 

(ri  (x,  t )rj(x,  f))V2  (r2(x,  f)r$(x,  f))1/2 

=  exp[— Af2/2r2j  (4.58) 

with 

At  =  B/2V 

Equation  (4.58)  expresses  the  ratio  between  the  INSAR  ocean  surface  cross-covariance 
of  complex  reflectivity  and  the  two  individual  SAR  covariance  of  complex  reflectiv¬ 
ities. 

Coherent  radars  such  as  SAR  are  inherently  subject  to  speckle  or  fading 
noise.  By  averaging  neighboring  pixels  the  speckle  noise  is  reduced  at  the  expense 
of  nominal  resolution.  This  operation  is  similar  to  the  multiple  look-operation,  in 
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which  the  processor  optimal  (maximum)  integration  time  is  divided  into  shorter 
integration  periods,  and  by  averaging  similar  images  with  a  time  shift,  the  fading 
noise  is  reduced  at  the  expense  of  resolution. 

In  the  data  analysis,  relatively  large  identical  areas  (>  103  pixels)  for  the 
two  SAR  and  INSAR  intensity  images  were  chosen  for  averaging.  The  noise  power 
is  reduced  by: 


where  N  is  the  number  of  averaged  pixels.  The  speckle  and  system  noise  are  further 
reduced  by  normalizing  the  ratio  K  in  Equation  (4.58)  by  a  similar  ratio  K0  obtained 
from  a  stationary  scene,  like  land,  which  fortunately  was  in  the  images.  The  scene 
coherency  of  land  varies  between  fields  with  objects  like  plants  and  trees  moving 
with  the  wind,  and  urban  regions,  which  should  be  absolutely  stationary  and  hence 
should  ideally  have  infinite  scene  coherence  time.  Due  to  speckle  noise  and  system 
errors  the  apparent  coherency  of  urban  regions  is  not  infinite;  thus  normalizing 
I<  from  Equation  (4.58)  with  K0  obtained  from  the  urban  regions  of  Marina,  will 
further  reduce  the  noise  resulting  in  a  better  estimation  of  rc  using  this  method. 
Based  on  this  assumption,  Equation  (4.58)  is  modified  to  be: 

K/K0  =  exp[-At2/2r2]  (4.60) 


and  the  scene  coherence  is  taken  to  be 


_ At 

^-21n(iiy/<o) 


(4.61) 


Equation  (4.61)  demonstrates  the  sensitivity  of  the  scene  coherence  time  estimation. 
A  change  in  the  ratio  K/Kq  of  2%  will  change  the  estimation  in  the  scene  coherence 
time  by  almost  20%  (see  Figure  4.4).  The  estimates  of  rc  using  the  above  experi¬ 
mental  method  yields  a  more  credible  L-band  scene  coherence  time  0  (100  msec), 
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V.  REFRACTION  SHOALING  AND 
DIRECTIONAL  SPECTRA  OF  NEARSHORE 

OCEAN  WAVES 

A.  INTRODUCTION 

This  chapter  presents  a  brief  review  of  linear  ocean  surface  waves  theory  in 
shallow  waters.  A  numerical  refraction  model  (Dobson,  1967),  followed  by  wave 
spectral  analysis  and  measurement  methods,  are  described.  The  behavior,  analysis, 
and  in  situ  observation  of  ocean  surface  waves  are  well  developed  subjects.  How¬ 
ever,  many  of  the  concepts  presented  are  difficult  to  interpret  in  terms  of  coherent 
and  noncoherent  radar  ocean  surface  imagery.  The  content  of  the  present  chap¬ 
ter  is  essential  for  interpretation  of  the  experimental  results  and  interdisciplinary 
understanding  of  the  issues  involved. 

B.  TRANSFORMATION  OF  NEARSHORE  OCEAN 

WAVES 

When  waves  propagate  from  deep  to  shallow  water  they  are  transformed  as 
they  begin  to  ‘feel’  the  bottom.  The  transformation  includes  shoaling,  refraction, 
diffraction,  reflection,  and  dissipation.  In  the  present  analysis  only  shoaling  and 
refraction  phenomena  are  considered,  since  we  ignore  reflection  (beach  slope  less 
than  1/20)  and  assume  energy  is  conserved  at  steady  state  and  is  not  translated 
lateraly  along  the  wave  crest  (ignoring  dissipation  and  difraction).  The  experimental 
site  at  Marina  provides  an  opportunity  to  study  wave  refraction  and  shoaling  as 
imaged  by  INSAR  and  ground  truthed  by  a  shallow  water  array  ofF  Marina. 
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1.  Wave  Refraction  and  Shoaling 


Starting  with  the  wave  energy  flux  balance  assuming  an  idealized  situa¬ 
tion  in  which  wave  energy  is  conserved  at  steady  state  and  not  transferred  lateraly, 
yields  the  simplified  energy  flux  equation: 

F 

(5.1) 


^k  +  ^  =  0 

dx  r  dy 


where  F (x,y)  is  the  wave  energy  flux  and  is  constant.  Given  a  monochromatic 
deep  water  gravity  wave  field  propagating  shoreward  as  depicted  in  Figure  5.1,  the 
following  energy  flux  balance  applies 

EoCgobQ  =  E\Cg\b\  (5.2) 


where  E0,E\  are  wave  energy  Cg0,Cgi  are  group  velocities  and  b0,  bx  are  distances 
between  wave  rays  in  deep  and  shallow  water  respectively.  We  further  discurs  the 
derivation  of  a  wave  height  relationship  for  shoaling  waves  from  Equation  (5.2)  based 
on  linear  wave  theory. 

2.  Refraction 

Assuming  the  gravity  waves  are  linear  and  dispersive  yields  the  relation¬ 
ship. 

u2  =  gkta.nh  kh  (5.3) 


where  u>  is  the  wave  frequency,  g  is  the  acceleration  of  gravity,  k  =  2i r/L  is  the 
wavenumber  magnitude  ( L  is  the  wavelength)  and  h  is  the  water  depth.  The  group 


velocity  for  a  surface  gravity  wave  train  is  defined  as 

duj 


Cn  = 


Ok 


(5.4) 


Substituting  Equation  (5.4)  into  Equation  (5.3)  yields  the  group  wave  velocity  for 


gravity  waves 


r  _  c  A  , 

9  2  \  sinh  2 kh 


(5.5) 
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Wove  roys 


Figure  5.1:  Wave  refraction  and  shoaling, 
where  C  is  the  wave  phase  velocity. 

The  instantaneous  surface  elevation  of  a  monochromatic  wave  is  given 


by: 


t/(x,  t)  =  asin(k  •  x  —  tat)  (5.6) 

where  a  is  the  wave  amplitude  and  x  =  (a :,y)  is  the  position  vector.  The  variance 
of  this  wave  is 

=  j  jQ  ^?(x,  0  -  dx=J  (5-7) 

where  L  is  the  ocean  wavelength.  Therefore  the  energy  of  a  monochromatic  wave 
in  an  irrotational  and  incompressible  fluid  given  by: 

E  =  \  P9a2  =  ^  P9U7  (5-8) 

where  H  —  2a  is  the  wave  height. 
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Substituting  Equation  (5.8)  into  Equation  (5.2)  yields  the  wave  height 
variation  equation  for  shoaling  waves  (Mei,  1980) 

■p  1 1/2  r,  1/2 

H  =  H°C%  T  =Hok,kR  (5.9) 


where  k,  is  the  shoaling  coefficient 


(5.10) 


and  kn  is  the  refraction  coefficient 


[h  11/2 

kn=[j\  •  (5-11) 

Within  the  limits  of  linear  wave  theory,  the  wave  frequency  must  re¬ 
main  constant  as  waves  propagate  from  deep  water  ashore.  Using  this  approximation 
in  the  dispersion  relationship  rewritten  in  a  wavelength  form  results  in 


L' i9h  (ir) = L° igh  (x1 


(5.12) 


Equation  (5.12)  illustrates  the  variation  in  wave  length  with  depth  for  shoaling  waves 
to  preserve  the  invariant  wave  frequency.  For  deep  water  waves 


(5.13) 


and  Equation  (5.13)  can  be  rewritten  as  (Shuchman  and  Kasischke,  1981) 


L\  =  L0  tgh 


(5.14) 


For  straight  and  parallel  of  bathymetric  contours,  a  simplified  refraction 
formula  based  on  Snell’s  low  can  be  stated  as: 


sinai 

Oi  =  Oo  — - 

sin  qq 


(5.15) 
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where  C  =  w/k  is  the  wave  phase  speed.  Using  Equation  (5.3)  to  compute  the  phase 
speed  and  substituting  the  result  in  Equation  (5.15)  results  in  the  refraction  angle 
a  between  the  wave  rays  and  bottom  contours  for  waves  propagating  from  deep  to 


shallow  water 


aj  =  sin-1  sin  c*o  tan  h 


(¥)] 


(5.16) 


In  the  data  analysis  described  in  Chapter  VII,  the  shoaling  effect  is 
estimated  by  using  Equations  (5.12)  and  (5.14),  and  the  refraction  effect  is  estimated 
using  Equation  (5.16)  for  straight  and  parallel  contours  only.  An  attempt  to  estimate 
refraction  and  shoaling  of  coastal  gravity  waves  using  SEASAT  SAR  data  was  done 
by  Shuchman  and  Kasischke  (1975).  The  method  presented  in  this  study  is  seen  to 
be  more  accurate. 

3.  The  Ray  Trace  Model 

For  complex  bathymetric  contours  which  are  not  straight  and  parallel  to 
the  shoreline,  a  linear  numerical  refraction  and  shoaling  model  by  Dobson  (1967)  is 
used  in  our  analysis.  Following  the  approach  established  by  Munk  (1952)  on  wave 
refraction,  Dobson  developed  a  linear  refraction  model  using  the  characteristics  of 
the  propagation  of  curving  wave  rays.  Considering  the  geometry  illustrated  in  Figure 
5.2  the  ray  curvature  equation  in  Cartesian  coordinates  is  given  by: 


n/  1  [•  dCix>y)  ,>sdC{.x,y) 

R(x'V'l)=C&y)  [SmQ(i)“s; - cosa<!>-0T“ 

and  the  ray  separation  factor  /?  is  stated  by  Dobson  to  be: 

+«(W=0 


(5.17) 


(5.18) 


dC(x,y) 

p(x,  y,  <)  =  2  cos  a — ^ - bsina 


dC(x,y ) 
dy 


(5.19) 
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Figure  5.2:  Refraction  geometry  and  notation  for  Dobson’s  model. 


,  J.2  d*C(x,y)  n.  d2C(x,y)  ,  2  0*C(*,y) 

q(x,y,t)  =  C(x,y)  sura — q„2  ~  2sinQr cos  a  - h  cos  a — — 

(5.20) 

The  wave  ray  path  defined  by  R  and  a  is  computed  with  Equations  ( 17)— (20), 
provided  the  wave  phase  speed  C(x,y)  is  known.  C(x,y)  is  constructed  grid  to  grid 
over  the  ocean  bathymetry  as  propagating  to  the  required  position.  For  more  details 
of  this  method  see  reference  (Dobson,  1967). 

In  this  study,  the  bathymetry  used  in  Dobson’s  refraction  model  is  based 
on  data  collected  originally  by  NOAA  and  projected  onto  a  six  second  (lat,  long) 
Modified  Universal  Transverse  Mercator  (UTM)  grid.  The  data  was  screened  and 
interpolated  using  linear  triangulation  to  fill  missing  points  onto  a  200  meter  rect¬ 


angular  grid  (Kalid,  1989). 
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C.  OBSERVATION  METHODS  OF  OCEAN  WAVE 
SPECTRA 

The  relation  between  the  distribution  of  energy  in  sea  surface  waves  and  the 
direction  of  propagation  is  of  great  interest  from  the  forecasting  point  of  view  and 
the  understanding  the  generation  processes.  For  the  case  presented  in  this  study, 
shoaling  wave  transformations  such  as  diffraction  ,  refraction  and  reflection  (see 
next  section)  are  greatly  influenced  by  the  directional  spectral  characteristics  of  the 
waves.  However,  measurements  of  directional  wave  spectra  requires  greater  effort 
than  that  of  frequency  power  spectra.  While  the  latter  can  be  obtained  from  a  wave 
record  at  a  single  point,  the  former  requires  simultaneous  recording  of  several  wave 
components. 

The  measuring  techniques  which  are  known  to  be  applied  can  be  classified 
into  two  categories  (Goda,  1985) 

1)  In  situ  measurement  method 

•  Wave  gauge  array 

•  Directional  buoy 

•  Two  axis  current  meter 

2)  Remote  sensing  method 

•  Optical  technique  (sterophotogrametry,  holography) 

•  Microwave  technique 

In  Chapter  VII,  frequency  and  directional  spectra  obtained  from  a  wave  gauge 
array  are  compared  with  two-dimensional  wavenumber  spectra  obtained  by  a  remote 
sensing  interferometric  SAR.  Therefore,  in  the  present  section  a  synopsis  of  the  basic 
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theoretical  and  observational  background  of  direct  and  remote  sensing  methods  are 
described. 

1.  Ocean  Surface  Wave  Field  Spectra 

The  characteristics  of  wave  field  properties  can  be  determined  from  the 
observed  sea  surface  elevation  and  the  irrotational  surface  velocity.  From  linear  wave 
theory  there  is  a  direct  relationship  between  the  sea  surface  elevation,  velocity  and 
pressure  spectrum.  The  horizontal  component  of  the  wave  orbital  velocity  can  be 
described  by  a  product  of  the  ocean  surface  elevation  and  a  transfer  function 


// „(w): 

u{t)  =  Hu{u)  r/(x,<) 

(5.21) 

where 

?;(x,  i)  =  a  cos(k  •  x  —  wi) 

(5.22) 

and 

„  ,  ,  cosh  k(h  +  z) 

=  “  sinh  kh  ' 

(5.23) 

where  z 

is  the  depth  of  observation.  The  wave  amplitude  a 

is  assumed  to  be  small 

compared  with  wavelength  and  depth  of  water  as  required  by  linear  theory.  Similarly, 

the  vertical  velocity  component  is  given  by 

w(<)  =  H*(u>)  r)(x,t) 

(5.24) 

where 

TT  ,  s  .  sinh  k(h  +  z) 

“  sinh  kh 

(5.25) 

The  pressure  wave  form  at  depth  z  is  taken  to  be 

P{t)  =  Hp{u) 

(5.26) 

where 

T.  .  .  cosh  k(h  +  z) 

,Ip{ul^PS  cosh  kh 

(5.27) 
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In  the  direct  method,  time  series  records  are  acquired  from  the  random 
ocean  processes.  Assuming  a  stationary  ocean  surface,  the  classical  analysis  ap¬ 
proach  is  the  spectral  technique,  i.e.,  linear  superposition  of  an  infinite  number  of 
sinusoidal  components  with  small  amplitude  and  slope.  Considering  the  water  col¬ 
umn  as  a  constant  parameter  linear  system,  the  input-output  frequency  is  invariant. 
For  input  surface  elevation  spectra,  Equations  (5.21),  (5.24),  and  (5.26)  can  be  ex¬ 
pressed  in  terms  of  input  output  power  spectra  related  by  the  corresponding  linear 
frequency  invariant  transfer  function  H(u>): 

S„,,M  =  |77(w)|  (5.28) 

Equation  (5.28)  reflects  the  fundamental  idea  of  measuring  the  ocean  wave  frequency 
power  spectra.  The  random  sea  surface  elevation  is  considered  as  the  input  and  the 
pressure  or  velocities  as  output.  Next  we  discuss  relevant  statistical  characteristics 
of  the  sea  surface  elevation  related  to  the  ocean  wave  spectra. 

The  random  ocean  surface  is  characterized  by  the  autocovariance,  p( r,  r), 
of  the  surface  displacement,  77  as: 

p(r,  t)  =  (*?(*>  0  >?(x,  r  +  r,  t  +  r))  (5.29) 

where  r  is  the  displacement  vector. 

The  ensemble  average  of  the  autocovariance  surface  elevation  in  Equa¬ 
tion  (5.29)  defines  the  decorrelation  time  and  space  scales  of  the  sea  surface.  The 
Wiener-Khinchin  theorem  states  that  for  any  ergodic  stationary  random  process 
(e.g.,  ? ?(x,r))  with  zero  mean,  the  wavenumber-frequency  spectrum  is  given  by 
the  Fourier  transform  of  the  autocovariance  function  p(r,  r)  of  the  random  process 
t?(x,  t). 

S(k,u>)  =  -^p  J  J  P(r,r)  exp[-i(k-r-u>7)]dr  dr  (5.30) 
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with  the  inverse  Fourier  transform 


/>(r.r)  =  /  /  £(k  ,u>)  exp[i(k- r -u>r)]<ik  du  (5.31) 

The  mean  square  of  surface  elevation  (J72(r,  r))  is  obtained  from  Equation  (5.30)  for 
r  =  0  and  r  =  0.  This  variance  term  is  proportional  to  the  wave  energy. 

Under  the  above  approximations  the  Fourier  transform  of  the  autoco¬ 
variance  of  the  surface  elevation  is  computed  as  a  function  of  wavenumber  and 
frequency,  with  position  and  time  as  parameters. 

2.  Relation  to  Measurement  Methods 

The  computation  of  the  frequency  wavenumber  spectra  requires  simul¬ 
taneous  measurements  of  surface  displacement  in  space  and  time.  This  requires  an 
array  of  wave  sensors  and  an  enormous  data  processing  effort.  The  omni-directional 
frequency  spectrum  is  obtained  with  less  effort  and  therefore  more  common,  since 
it  requires  only  a  single  wave  sensor  at  a  point,  collecting  temporal  data  with  no 
spatial  information  about  the  wave  field.  At  a  given  position  r  =  0,  the  frequency 
spectra  is  given  by 

S{u)  =  J  5(k,u>)  dk  —  —■  J  p{0,r)  exp[zW]c/r  (5.32) 

Since  only  positive  circular  frequencies  have  physical  meaning  and  S(u>)  is  symmetric 
about  w  =  0,  Equation  (5.32)  can  be  rewritten  as: 

1  r+°° 

S(u)  =  —  /  p(r)  cos  ut  dr  (5.33) 

7 r  Jo 

Since  the  total  energy  for  a  given  sea  state  is  unique,  the  integration  areas 
of  the  frequency  and  wavenumber  spectra  are  related  by  the  Jacobian,  yielding  the 
unidirectional  wavenumber  spectrum 

= 


(5.34) 


where  k  —  |k|  and  du>/dk  =  Cg  the  group  wave  velocity.  Using  the  linear  dispersion 
relation  u  =  u(k)  in  computing  Cg,S(k)  is  obtained.  This  transformation  can  not 
be  used  in  the  high  wavenumber,  short  waves  part  of  the  spectrum,  since  the  intrinsic 
frequency  of  the  shortest  waves  w0(&)  are  Doppler  shifted  due  to  modulation  by  the 
underlying  dominant  waves  and/or  currents  on  which  the  short  waves  are  borne. 
Therefore  the  circular  frequency  of  the  short  waves  is  given  by 

uj(k)  =  u0(k)  +  k  •  (U0b(x,  t)  +  UCb(x, 0)  (5.35) 

where  UoB(x,t)  is  the  slowly  varying  orbital  velocity  and  Ucn(x,<)  could  be  any 
underlying  varying  currents.  Equation  (5.35)  illustrates  the  nonuniqueness  of  the 
relationship  between  frequency  and  wavenumber. 

For  high  resolution  microwave  radar  sensors  the  backscattering  energy 
is  determined  by  the  short  scale  (cm)  surface  roughness  which  is  borne  on  an  ad- 
vecting  surface.  Equation  (5.35)  replaces  the  linear  dispersion  relationship  between 
frequency  and  wavenumber  (Equation  (5.3)). 

3.  Wave  Directional  Spectrum 

As  noted  in  the  literature  review  (Equation  (2.7))  the  normalized  radar 
cross-section  a0  is  directly  proportional  to  the  two  dimensional  wavenumber  spectrum. 

t/>( k)  =  J  S{k,w)dw  (5.36) 

The  instantaneous  autocovariance  can  be  expressed  as: 

/?(r,0)  =  J  J  ip(k  exp[-i(k  •  r)]dk  (5.37) 

and  inversely 

V»(k)  =  2  / pW  exp[{(k  •  r)}dr  (5-38) 
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To  obtain  ^>( k ),  one  needs  spatial  observations  of  the  surface  elevation  autocovari¬ 
ance  (Equation  (5.38)).  Attempts  to  measure  this  phenomena  are  still  under  exten¬ 
sive  research  using  remote  sensing  techniques  such  as  sterophotogrametry  (Shemdin, 
et  al.,  1988)  laser  optics  (Keller  and  Gotwols,  1983)  or  Doppler  radars  (Plant  and 
Keller,  1983). 

The  two-dimensional  wavenumber  spectrum  is  also  claimed  to  be  com¬ 
puted  from  SAR  images.  However,  there  is  still  doubt  about  the  precise  mechanisms 
whereby  waves  are  imaged  by  SAR.  Therefore,  the  question  remains  whether  the 
SAR  image  intensity  spectrum  bears  much  resemblance  to  the  directional  wave  spec¬ 
trum.  In  this  study,  a  new  approach  is  suggested  for  obtaining  a  two-dimensional 
wavenumber  spectrum  using  interferometric  SAR  imagery.  This  method  has  a 
clearly  observed  imaging  mechanism  based  on  the  measurement  of  the  varying  or¬ 
bital  velocity  Uob{x,1)  of  the  dominant  surface  ocean  waves.  The  two-dimensional 
wave  number  spectrum  of  INSAR’s  phase  image  represents  the  spatial  distribution 
of  the  surface  velocity  power  in  the  viewing  direction.  In  describing  this  spectrum, 
assume  for  the  moment  a  monochromatic  dominant  wave  oriented  at  an  arbitrary 
angle  relative  to  the  flight  direction  (see  Figure  4.1)  and  an  incidence  viewing 
angle  7.  The  surface  elevation  of  this  wave  as  seen  by  the  imaging  platform  is  given 
by: 

t?(x, t)  =  a  g{'i,(}>)  cos(k-x  —  ut)  (5.39) 

where  a  is  the  wave  amplitude  and 

g{  7,  <f> )  =  (cos2  7  +  sin2  7  sin2  (5.40) 

After  removing  the  current  and  phase  velocities,  the  residual  imaged 
radial  velocity  component  UoBr  can  be  computed.  Using  deep  water  linear  wave 
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theory 


tW(x)  =  ~gf'  =  auJ  9{r,  $)  sin(k  -x-ut) 


(5.41) 


Then  the  two-dimensional  wavenumber  spectra  for  the  orbital  velocity  is 

i’UoBA^K)  ~  p  /  /  ,  V0Br(x,y)  exp[-i(kxx  +  kvy)}  dxdy  (5.42) 

XT  J-L/2  J—L/2 

The  spectrum  described  in  Equation  (5.42)  is  tested  against  directional  spectra 
obtained  from  a  wave  pressure  array.  A  synopsis  review  on  the  theory  of  in  situ 
measurement  of  directional  spectra  is  described  next. 

Longuet-Higgins,  Cartwright  and  Smith  (1963)  developed  a  directional 
wave  bouy  based  on  measuring  the  angles  of  pitch  and  roll  and  the  heaving  accel¬ 
eration  of  the  ocean  surface.  The  heave  is  directly  related  to  the  surface  elevation, 
while  the  pitch  and  roll  are  related  to  the  two  horizontal  orthogonal  components  of 
wave  surface  slope  namely,  dy/dx ,  d~}/dy.  The  directional  spectrum  S(f,0)  can  be 
expressed  as  a  product  of  the  one-dimensional  frequency  power  spectrum  S(f)  and 
an  angular  spreading  function  D(/,  a) 

S(/,a)  =  S(/)£(/,a)  (5.43) 

where  a  is  the  wave  direction.  D(f,  a)  is  associated  with  the  distribution  of  energy 
in  the  frequency  power  spectrum,  over  direction  and  has  a  unit  area  over  direction, 

f  D(f,a)  da  =  1  (5.44) 

Jo 

An  objective  of  measurements  is  to  determine  the  characteristics  and  criteria  in 
selecting  D(f,a )  which  most  closely  match  the  observed  angular  distribution  of 
wave  energy  density.  Using  the  most  common  in  situ  measurement  technique  -  a 
floating  buoy  or  slope  array  of  four  sensors  in  a  square,  one  has  to  address  their 
limited  spatial  capabilities  in  defining  directionality,  namely  D(f,a). 
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Following  Longuet-Iliggins  et  al.  (1963),  the  directional  spreading  func¬ 
tion  £>(/,  a)  can  be  expanded  as  a  Fourier  series: 

/)(/, a)  =  ^  Q  +  J3[a„ cos(nan)  +  bns'm(nan)]j  (5.45) 

For  the  previously  discussed  bouys  or  slope  arrays  which  measure  only  the  surface 
slope  and  displacement,  only  the  first  two  harmonics  can  be  estimated.  The  hori¬ 
zontal  sea  surface  slopes  are  obtained  by  taking  the  corresponding  spatial  derivative 
of  the  sea  surface  elevation  (Equation  5.6) 

Vx{t)  =  =  (—ik  cos  or)  »/(<)  =  Hx(f)  i](t)  (5.46) 

Vy(l)  =  =  (~ik  cos  or)  =  Hv{f)  rj(t)  (5.47) 

Five  Fourier  coefficients  can  be  determined  from  the  co(C)  and  quad  ( Q )  spectra  of 
the  surface  elevation  and  surface  slope.  The  spectrum  of  tjx(t)  is  given  as 

^vxvx  (/)=  I**  \IUf)\2S(f,a)da  (5.48) 

Jo 

Substituting  Equations  (5.43),  (5.46)  into  Equation  (5.48)  yields 

SVir,M)  ~  f  k2  cos2  aS{f)D(f,a)  da  =  k2S{f)  f  cos2  aD(f,a)  da  (5.49) 
Jo  Jo 


Similarly  the  rest  of  cross-spectral  quantities  are  obtained 

Cw,(f)  =  fc2S(/)  J*’  sin2  a  D(f,  a)  da  (5.50) 

£,„(/)  =  S(f)  r  D(f,a)da  =  S(f)  (5.51 ) 

Jo 

c„„(/)  =  USU)  jf  cos  a  sin  aD(f,a)  da  (5.52) 

Qmx(f)  =  kS(f)  f2*  cos  a  D{f,a)  da  (5.53) 

Jo 

Qnv„(f)  =  kS{f)  [  sin  a  £>(/,«)  da  (5.54) 

J  0 
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From  Equations  (5.49)-(5.54)  Longuet-Higgins  et  a!.  (1963)  computed  the  Fourier 


coefficients  of  £>(/,  a)  as  follows: 

a0  =  —  [  D(f,a)  da—-  (5.55) 

7T  JO  7T 

ai  =  -[  D(/,a)  cos  a  da  =  Qvv.if)  (5-56) 

7T  JO  ,  7TK 

b‘  =  ~J0  DU,a)  sin  a  da  -  ~  Qm,(!)  (5.57) 

*>  =  \ J*  W.c)  cos 2a  da  =  ^  (C*„(/)  -  C„,„(/))  (5.58) 

5a  =  i  D(f,a)  sin  2a  <ia  =  -|j  C„,„,(/)  (5.59) 


Longuet-Higgins  et  al.  (1963)  used  these  five  Fourier  coefficients  to  de¬ 
scribe  an  assumed  cosine  bell  angular  energy  distribution  function.  This  method 
can  only  describe  a  unimodal  energy  distribution  with  limited  resolution,  which  is 
frequently  not  realistic.  These  coefficients  can  be  represented  by  an  infinite  number 
of  angular  spreading  functions.  It  is  therefore  desirable  to  seek  a  technique  (if  pos¬ 
sible)  to  match  the  coefficients  obtained  from  the  data  to  a  more  realistic  angular 
distribution  function. 

Grauzinis  (1990)  improved  the  Longuet-Higgins  et  al.  (1963)  method 
in  matching  analytically  a  more  realistic  angular  distribution  function  D(a )  to  the 
obtained  Fourier  coefficients.  The  technique  is  based  on  representing  the  Fourier 
coefficients  of  D(a)  as  the  superposition  of  a  weighted  sum  of  symmetric  unit  area 
density  functions  oriented  at  different  direction 

a*  +  ibk  =  Ew;  mk{j)  cxp(jfc/?;)  (5.60) 

3 

where  w;-  is  the  weighting  parameter,  m*  is  the  canonical  coefficient  of  the  even  unity 
spreading  function  and  fy  is  the  angular  orientation  of  component  j.  For  a  typical 
pressure  array  (four  sensors)  or  pitch  and  roll  buoy,  the  obtained  cross  spectral  wave 
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quantities  result  in  the  Fourier  coefficients  of  Equation  (5.60)  through  the  second 
order  and  are  taken  to  be 


+  ib i  =  mi  [wi  exp(t/?i)  -f  w2  exp(i/?2)] 

(5.61) 

a2  + 1&2  =  m2  [w!  exp(i2/?i)  4-  w2  exp(z'2/?2)j 

(5.62) 

To  solve  the  system  represented  by  Equations  (5.61),  (5.62),  further  as¬ 
sumptions  have  to  be  made.  The  most  sufficient  assumption  made  by  Grauzinis 
(1990)  is  that  the  two  canonical  coefficients  are  related  by  a  power  law  r 

m2  =  mj  (5.63) 

where  r  defines  the  broadening  (r  increase)  or  the  sharpness  (r  decrease)  of  the 
distribution  function  defined  initially  by  the  canonical  coefficients.  For  a  unimodal 
wave  field  and  m  =  0.9  and  r  =  4,  the  half  power  distribution  is  12  degrees,  which 
is  about  one  third  of  the  Longuet-Higgins  et  al.  (1963)  cosine  bell  resolution.  The 
improvement  of  this  technique  relative  to  the  Longuet-Higgins  et  al.  (1963)  method 
is  the  variability  in  the  sharpness  of  the  model  distribution  function,  which  enables 
it  to  be  matched  with  better  directionality  to  the  observed  energy  distribution  func¬ 
tion.  Also,  Equations  (5.61),  (5.62)  introduce  two  distributions  of  wave  energy  with 
different  orientations  and  weights.  This  allows  the  model  to  represent  a  bimodal 
wave  system  which  is  not  possible  by  the  Longuet-Higgins  et  al.  (1963)  method. 

This  exact  method  of  matching  the  Fourier  coefficients  to  a  wave  spread¬ 
ing  function  is  not  unique.  The  question  that  remains  is  how  to  optimize  the  selec¬ 
tion  of  the  distribution  that  best  matches  the  real  data.  Cartwright  (1963)  suggests 
that  a  criteria  to  improve  the  selection  of  the  distribution  function  can  be  obtained 
from  the  ratio  of  the  magnitude  of  the  Fourier  coeffients.  Grauzinis  (1990)  further 
suggests  a  beam  separation  index  criteria  to  indicate  the  best  distribution  function 
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fitting.  This  criteria  aids  the  selection  of  the  canonical  coefficients  which  fulfills  the 
analytic  solution  to  give  zero  beam  separation  from  the  required  distribution.  For 
different  r’s,  the  bimodal  distribution  with  the  minimum  peak  separation  is  selected 
as  the  best  fit.  Despite  the  analytical  efforts  to  improve  the  directional  resolving 
wave  power  obtained  from  pitch  and  roll  buoys  or  pressure/current  limited  arrays, 
an  rms  resolution  of  4-6  degrees  appears  to  be  the  best  that  can  be  achieved  using 
this  method. 
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VI.  THE  MONTEREY  BAY  WAVE  FIELD 
AND  CURRENTS  ‘89  REMOTE  SENSING 

EXPERIMENT 

A.  INTRODUCTION 

The  primary  objective  of  the  experimental  part  of  this  study  is  to  verify  the 
feasibility  of  an  airborne  INSAR  method  to  reliably  spatially  sense  ocean  surface 
wave  and  current  fields.  The  verification  is  done  by  comparing  the  results  of  this 
unique  technique  with  simultaneous  tn  situ  conventional  temporal  and  local  spatial 
measurements  including  directional  wave  spectra  from  an  array,  local  wind  and  sur¬ 
face  drifters  measurements.  To  simplify  the  interpretations,  the  time  of  the  year 
and  the  site  of  the  experiment  were  selected  so  as  to  provide  near  linear  condi¬ 
tions  of  wave  imaging  (narrow  band  wave  spectrum  of  low  to  moderate  variance) 
and  a  moderate  sloping  beach  (approximately  20  m/km)  allowing  sufficient  locally 
homogeneous  conditions  for  linear  spectral  analysis. 

The  relatively  shallow  water  provides  an  opportunity  to  study  for  the  first 
time  wave  refraction  and  shoaling  using  INSAR.  Four  partially  overlaped  areas  (~ 
12  km  x  6  km)  were  imaged  by  INSAR  from  four  different  directions  centered  on 
the  shallow  water  wave  array  ofT  Marina.  Results  from  three  legs  are  reported  in 
this  study. 

Overall,  three  field  tests  were  performed  in  the  vicinity  of  Monterey  Bay: 

•  On  June  1,  1989,  several  prototypes  of  passive  radar  drifters  were  'ted  at 
sea  to  measure  ocean  surface  currents  using  Lagrangian  bouys.  The  drifters 
were  designed  and  built  at  the  Naval  PostGraduate  School  and  were  tested 
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from  the  Research  Vessel  (R/V)  Ricketts  in  the  vicinity  of  Moss  Landing  (see 
Appendix  C). 

•  On  July  23-24,  1989,  a  preliminary  coordinated  simultaneous  airborne  and 
ground  truthing  experiment  was  performed.  The  experimental  scene  was  cen¬ 
tered  on  the  wave  and  wind  stress  NDBC  pitch  and  roll  buoy  off  Monterey  Bay 
(see  Figure  6.2).  NASA/JPL  DC-8  radar  overflights  were  conducted  on  July 
24, 1989  between  1000  and  1040  Pacific  Day  Time  (PDT).  The  R/V  Point  Sur 
on  a  cruise  from  232200  to  242400  July  1989  acquired  ground  truthing  data 
for  comparison  with  the  following  near  ocean  surface  phenomena: 

-  Vertical  profile  of  current  vectors  within  the  mixed  layer  measured  by 
Acoustic  Doppler  Current  Profile  (ADCP). 

-  Surface  currents  sensed  by  Lagrangian  floating  buoys  which  were  posi¬ 
tioned  by  Loran-C  and  Mini-Ranger  navigation  systems  (see  Section  E4). 

-  Wind  speed  and  direction  which  generates  the  short  ocean  surface  wave 
field  patterns  responsible  for  the  resonant  radar  backscatter  (Bragg  waves). 

—  Vertical  density  gradients  at  the  depth  between  10  and  250  m  measured 
by  Conductivity  Temperature  Depth  (CTD),  which  might  affect  surface 
dynamical  processes. 

Wind  stress  and  dominant  wave  time  series  records  were  acquired  from  the  NOAA 
pitch  and  roll  buoy.  Unfortunately,  the  Interferromatric  SAR  images  could  not  be 
obtained  from  this  experiment.  An  error  in  setting  up  the  fore  and  aft  radar  anten¬ 
nas  with  different  polarizations  (vertical  and  horizontal),  instead  of  being  identical 
(both  vertical  or  horizontal),  prevented  fringes  formation  essential  for  creating  an 
interferogram. 
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Based  on  these  preliminary  experiments,  the  experimental  plan  was  modified 
for  a  following  experiment.  The  field  test  site  was  moved  to  an  area  close  to  the  shore 
to  include  stationary  land  in  the  INSAR  imaging  scene  to  assist  aircraft  attitude 
corrections  in  addition  to  the  aircraft  Inertial  Navigation  System  (INS)  (see  Chapter 
III).  A  second  modification  was  to  replace  the  passive  radar  drifters  with  thin  floating 
plywood  sheets  (see  section  F3)  to  minimize  windage  effects.  The  following  and  final 
experiment  for  this  study  was  carried  out  on  September  8,  1989  in  Monterey  Bay 
centered  on  the  shallow  water  wave  array  off  Marina  Beach.  An  illustration  of  the 
July  and  September  experimental  scenarios  is  depicted  in  Figure  6.1.  A  detailed 
description  of  the  Marina  experimental  site  is  given  in  section  B,  followed  by  a 
description  of  the  Marina  experiment  in  section  C.  The  environmental  conditions 
for  the  experiment  are  given  in  section  D. 

B.  THE  EXPERIMENTAL  SITE 

The  Monterey  Bay  is  characterized  by  variable  bathymetric  features  including 
the  largest  submarine  canyon  in  the  western  hemisphere  bisecting  the  Bay  (Figure 
6.2).  These  bottom  features  cause  substantial  perturbations  on  the  near  shore  wave 
fields  and  tidal  currents.  The  experimental  site  for  data  acquired  during  this  study 
is  in  the  nearshore  region  of  Marina  within  Monterey  Bay.  This  region  is  shadowed 
by  Point  Santa  Cruz  in  the  northwest  and  Point  Pinos  in  the  southwest. 

The  Marina  site  was  selected  because  of  the  location  of  a  wave  array  for 
ground  truthing  and  the  required  inclusions  of  land  in  the  scene.  Also,  the  site 
has  a  relatively  continuous,  though  moderate,  alongshore  bottom  slope  (about  20 
m/km)  resulting  in  a  relatively  small  but  sufficient  region  of  locally  homogeneous 
bathymetry.  This  fact  constrained  the  spectral  analysis  to  be  applied  to  small 
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Figure  6.1:  An  artistic  illustration  of  experiment  scenario  conducted  on 
July  24  and  September  8,  1989. 


spatial  regions  to  meet  the  required  homogeneous  conditions.  This  area  presents  an 
interesting  problem  in  predicting  directional  wave  spectra. 

The  sloping  bottom  introduces  additional  various  dynamic  ocean  features 
like  gradients  in  the  tidal  flows.  The  present  study  concentrates  on  the  dominant 
surface  waves  spectra  and  currents.  Shoaling  waves  undergo  refraction,  reflection 
and  diffraction  influenced  by  the  wave  approach  angle  and  bathymetry.  The  complex 
bathymetry  of  Monterey  Bay  is  shown  in  Figure  6.2,  and  the  detailed  bathymetry 
near  Marina  in  Figure  6.3.  Wave  refraction  at  the  Marina  array  site  is  shown  in 
Figure  6.4  for  different  frequencies  (0.03  Hz  to  0.13  Hz).  Linear  refraction  theory 
(see  Chapter  V)  predicts  that  lower  frequency  waves  (longer  waves)  experience  more 
refraction.  All  possible  incoming  ray  traces  with  wave  angles  at  increments  of  0.2° 
using  Dobson’s  model  (see  Dobson,  1967)  are  plotted  from  shallow  to  deep  water. 
High  frequency  waves  of  at  least  0.11  Hz  are  limited  to  directions  from  234°-311°, 
whereas  lower  frequency  waves  may  arrive  aiso  from  southerly  directions  since  they 
experience  more  refraction. 

An  important  feature  in  the  experimental  site  is  the  Salinas  River,  which  was 
blocked  from  the  sea  at  this  time  of  the  year  (June-September)  providing  a  natural 
closed  basin.  This  allowed  the  study  of  Bragg  backscattering  waves  only,  without 
the  contamination  by  interactions  with  long  waves  or  any  underlying  currents. 

C.  THE  MARINA  EXPERIMENT 

The  Marina  experiment  consisted  of  four  SAR  and  INSAR  overflights  in  the 
interferometric  mode  on  September  8,  1989.  Ra.dar  data  were  acquired  about  an 
area  centered  on  the  shallow  water  wave  array  (see  Figure  6.5)  with  the  antenna 
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Figure  6.4:  Refraction  wave  ray  traces  off  Marina  Beach  at  different 
frequencies. 
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Figure  6.5:  Flight  pattern  and  scene  locations  for  radar  data  collection 
on  September  8,  1989  off  Marina  Beach. 
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beams  aimed  at  the  in  situ  wave  array  with  a  depression  angle  of  45°  and  mean  alti¬ 
tude  of  8250  m.  The  radar  overflights  were  conducted  around  1300  Pacific  Daylight 
Savings  Time  close  to  the  time  of  maximum  tidal  flood  current.  The  flight  pattern 
and  the  scene  locations  of  radar  data  collection  are  illustrated  in  Figure  6.5.  Some 
key  system  and  flight  characteristics  are  depicted  in  Table  6.1.  More  detailed  SAR 
system  characteristics  are  given  in  the  next  chapter. 

The  Lagrangian  floating  drifters  for  tracking  the  very  top  surface  currents 
were  deployed  and  tracked  by  the  R/V  Ricketts  with  the  aid  of  a  Piper  Cherokee 
Arrow  (PCA)  light  aircraft  (Figure  4.15).  A  portable  LORAN-C  enabled  the  light 
aircraft  to  locate  roughly  the  drifters  while  the  vessel,  and  hence  the  drifter  positions, 
were  recorded  using  a  high  accuracy  Motorola  Mini-Ranger  Falcon  484  system. 

The  actual  time  schedule  and  sequence  of  events  in  the  Marina  experiment 
on  September  8,  1989  were  as  follows: 

081000  Departure  of  R/V  Ricketts  from  Moss  Landing  home  base. 

081050  On  station  (Marina  shallow  water  array)  to  test  the  navigation  system 
(Mini-Ranger)  and  to  estimate  local  surface  currents  by  deploying  and 
positioning  two  drifters.  This  pre-measurement  was  done  in  order  to  op¬ 
timize  the  deployment  pattern  of  the  drifters  to  obtain  maximum  coverage 
of  variable  surface  currents  within  the  expected  airborne  radar  scene  for 
data  collection.  The  variability  of  the  bathymetry  was  part  of  the  concern 
in  ground  truthing  the  radar  surface  current  estimates  at  variable  depth. 
081120  Termination  of  initial  test,  vessel  steaming  westward  to  the  first  deep 
water  deployment  station. 

081145  -  081230  Buoy  deployment  from  west  to  east  in  three  distinct  groups, 
each  group  consisting  of  five  drifters  for  a  total  of  fifteen  (see  Figure 
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TABLE  6.1:  System  and  Flight  Data  for  the  Different  Aerial  DC-8  Air¬ 
craft  Coverage  During  Marina  Experiment,  September  8,  1989 


90° 

EAST¬ 

WARD 

LEG 

270° 

WEST¬ 

WARD 

LEG 

360° 

NORTH¬ 

WARD 

LEG 

Initial  data  collection 
time,  PDT 

12:40 

13:01 

12:51 

Aircraft  velocity,  m/s 

214 

216 

214 

Aircraft  altitude,  m 

8285 

8273 

8178 

Pulse  repetition  frequency  ,  Hz 

283 

285 

282 

Near  incidence  angle,  degrees 

17 

Nadir 

Nadir 

Far  incidence  angle,  degrees 

52.7 

50.2 

50.7 

Ground  pixel  range  spacing,  m 

8 

7 

6 

Azimuth  pixel  spacing,  m 

12.1 

12.1 

12.1 

Effective  coverage 
range  x  azimuth,  km2 

6  x  14 

5.1  x  14 

5.8  x  14 

6.6).  Group  one  drifters,  numbered  1  thru  5,  were  deployed  furthest 
ofTshore  approximately  3.4  nm  west  of  the  shore  where  the  mean  water 
depth  is  71  m.  This  group  was  positioned  four  times  during  the  three 
hour  observation  period.  Group  two  drifters,  6  thru  10,  were  deployed 
closer  to  the  shore,  approximately  1.6  nm  west  of  the  beach  where  the 
mean  water  depth  is  44  m.  These  were  positioned  three  times.  Group 
three  drifters,  11  thru  15,  were  deployed  less  than  0.6  nm  from  the  beach 
and  washed  ashore  before  a  second  position  could  be  obtained.  Within 
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each  group,  the  drifters  were  deployed  roughly  0.1  nm  apart  as  the  vessel 
steamed  eastward  toward  the  beach. 

081220  The  Piper  Cherokee  Arrow  light  aircraft  arrived  on  station  Monterey  Air¬ 
port  and  started  to  position  and  photograph  the  drifters  from  west  to 
east. 

081230  -  081515  Drifters  chased,  with  vessel  heading  west  to  position  the  drifters 
first  deployed.  The  initial  detection  of  the  Lagrangian  drifters  was  quite 
tedious.  Applying  Direction  and  Range  (DR)  estimates  based  on  the 
early  test  results,  the  drifters  were  not  detected.  The  light  aircraft  was 
very  helpful  in  detecting  the  drifters  initially,  although  the  radio  com¬ 
munications  broke  down.  Therefore  the  actual  communication  was  done 
by  turning  on  the  landing  lights  and  heading  towards  the  R/V  Ricketts 
once  the  Piper  overflew  a  group  of  drifters.  As  predicted,  the  onshore 
tidal  current  became  more  intense  close  to  the  time  of  midtide  range 
(1320),  yielding  a  stronger  eastward  component,  which  was  measured  by 
the  drifters.  The  second  DR  prediction  based  on  the  initial  physical  po¬ 
sitioning  was  close  to  the  actual  positions  of  the  drifters. 

081240  -  081310  SAR  overflight  acquired  data  in  the  INSAR  mode. 

081405  The  Piper  light  aircraft  heads  home  to  Monterey. 

081515  R/V  Ricketts  heads  home  to  Moss  Landing  after  retrieving  ten  of  the 
drifters  from  the  two  most  western  groups.  The  most  eastward  group  of 
drifters,  which  were  washed  ashore,  were  collected  four  days  later  from 
Marina  Beach. 
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Figure  6.6:  Drifter  deployment  positions  and  Mini-Ranger  reference  sta¬ 
tions  during  the  Marina  experiment,  September  8,  1989. 
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D.  ENVIRONMENTAL  CONDITIONS 

1.  Introduction 

The  environmental  conditions  which  occurred  during  the  Marina  exper¬ 
iment  on  September  8, 1989  were  typical  late  summer  conditions  for  Monterey  Bay. 
Light  westerly  winds  (less  than  4  knots),  surface  air  temperature  of  14-16°C  and 
clear  skies  up  to  1000  feet  blocked  by  a  stratus  cloud  deck,  were  observed  during 
the  1000-1500  Pacific  Daylight  Time  (PDT)  experiment  period  (see  Table  6.2). 

The  mild  local  meteorological  condith  ns  during  this  period  resulted  in 
almost  no  local  generation  of  sea.  The  sea  surface  was  predominantly  long  crested, 
narrow-band  waves  propagating  shoreward.  Therefore,  the  scene  was  linear  in  char¬ 
acter,  yielding  a  linear  Modulation  Transfer  Function  (MTF)  of  the  ocean  surface 
for  imaging  radars. 

Due  to  the  lack  of  predicted  tidal  current  information  for  Marina  Beach, 
the  time  of  maximum  tidal  current  was  estimated  using  predicted  information  of 
tides  at  Monterey  Harbor  and  Moss  Landing,  which  are  located  almost  in  equal 
distance  (17  km)  south  and  north  off  Marina  Beach.  Predicted  tidal  high  and  low 
water  times  (NOAA  1978  Tide  Tables)  of  Monterey  Harbor  and  Moss  Landing  are 
plotted  in  Figure  6.7.  Based  on  the  figure,  the  maximum  flood  tidal  current  was 
probably  at  1320  PDT. 

The  summer  climatological  atmosphere  surface  circulation  in  the  eastern 
North  Pacific  consists  of  a  high  east-west  pressure  ridge  stretching  from  west  coast 
California  to  almost  140°  E  (east  longitude).  The  ridge  extends  to  approximately  20° 
N  (north  latitude).  Below  20°  N,  the  atmosphere  has  tropical  characteristics.  This 
late  August  high  pressure  system  was  deformed  in  shape  and  orientation  (see  Fig¬ 
ure  6.8)  by  a  series  of  northwest  lows,  which  might  be  the  source  of  the  northwesterly 
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TABLE  6.2:  Environmental  Conditions  During  the  WAFICIN  Marina 
Experiment 


Time  and  date 

Air  temperature 
Surface  Water  temperature 
Relative  humidity 
Wind  direction 

Wind  speed 
Visibility 

Significant  wave  height 
Dominant  swell  periods 

Expected  time  of  maximum 
flooding  tidal  current 

Swell  directions 


1030-1500  PDT, 

September  8,  1989 

14  to  16°C 

~  13°C  (est.) 

~  60%  to  70%  (est.) 

Varying  between  northeast 
and  southeast 

<  2  m/s 

Up  to  1000  ft  7  miles. 
Above  1000  ft  overcast  by 
low  stratus  cloud  deck. 

<  0.6  m 

a.  9.1  sec 

b.  15.9  sec 
1320  PDT 

a.  East  Northeast 
(long  swell) 

b.  East  Southeast 
(short  swell) 
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Figure  6.7:  Monterey  Harbor  and  Moss  Landing  predicted  tides,  Septem¬ 
ber  7-9,  1989. 
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9.1  sec.  swell  with  travel  time  of  almost  4  days  observed  in  the  Marina  experiment 
on  September  8,  1989.  The  15.9  sec.  long  crested  westerly  swell  is  postulated  to  be 
generated  by  a  late  winter  southern  hemisphere  storm  located  at  40°  S  to  60°  S  in 
the  South  Pacific.  The  travel  time  for  the  swell  from  this  storm  to  Monterey  Bay  is 
almost  8  days  (see  Figures  G.8  and  6.10). 

2.  Evolution  of  Environmental  Conditions 

When  investigating  the  possible  sources  of  the  observed  swell  waves  in 
the  Marina  experiment,  three  different  sources  of  process  data  were  analyzed: 

•  Surface  winds  and  pressure  field  patterns  between  August  31  and  September 
8, 1989,  based  on  NORAPS/SAC  model  outputs  processed  by  Fleet  Numerical 
Ocean  Center  (FNOC)  in  Monterey  (see  Figures  6.8,  6.9). 

•  Wave  group  velocity  vectors  and  significant  wave  height  contours,  based  on 
Global  Spectral  Ocean  Wave  Model  (GSOWN)  outputs,  between  August  31 
and  September  8,  1989,  also  processed  by  FNOC  in  Monterey  (see  Figures 
6.10,  6.11). 

•  Geostationary  Operational  Environmental  Satellite  (GOES)  visible  imagary 
(Figure  6.12). 

3.  Surface  Winds  and  Pressure  Fields 

The  surface  pressure  field  pattern  on  August  31,  1989  in  the  eastern 
North  Pacific  was  in  good  agreement  with  climatology.  From  the  1st  of  September 
1989,  the  east-west  orientation  of  the  ridge  axis  shifted  to  a  northeast,  southwest 
orientation.  This  shift  was  caused  by  the  development  of  a  strong  low  pressure 
system  south  of  the  Aleutian  Islands  and  the  development  of  an  induced  trough 
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Figure  6.8:  Mean  sea  level  pressure  contours  and  wind  vectors  at  1700 
PDT  on  August  31,  1989. 
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Figure  6.9:  Mean  sea  level  pressure  contours  and  wind  vectors  at  1700 
PDT  on  September  8,  1989. 
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Figure  6.10:  GSOWN  wave  group  velocity  vectors  and  significant  wave 
height  contours  at  1700  PDT  on  August  31,  1989. 
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Figure  6.11:  GSOWN  wave  grou 
height  contours  at  1700  PDT  on  S 
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Figure  6.12:  GOES  visible  imagery  at  1331  PDT  on  September  8,  1989. 
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over  California  and  the  Pacific  Northwest.  These  pressure  field  conditions  were 
maintained  through  September  8,  1989  yielding  stable  conditions  with  light  winds 
parallel  to  the  coast  and  westerly  light  winds  between  35°  N  to  38°  N  (Figure  6.9). 
These  computed  wind  fields  based  on  NORAPS/SAC  model  arc  in  agreement  with 
the  observed  wind  tendancy  during  the  Marina  experiment. 

The  migration  of  the  pressure  ridge  northward  inducing  stable  light  wind 
conditions  over  the  west  coast  enabled  wave  systems  produced  by  northwest  Pacific 
storms  to  propagate  eastward  without  local  air-sea  disturbances.  As  mentioned 
previously,  these  developing  storms  appear  to  be  the  source  of  the  shorter  period 
swell  (9.1  sec.)  observed  on  September  8  in  the  Marina  experiment. 

4.  Satellite  Data 

Analysis  of  GOES  imagery  confirms  the  presence  of  an  induced  trough 
over  California  and  Pacific  Northwest  (see  Figure  4.12).  The  low  pressure  system 
south  of  the  Aleutian  Islands,  which  is  believed  to  be  the  source  of  westerly  winds, 
caused  the  1000  foot  ocean  stratus  deck  to  approach  the  central  and  North  California 
coast.  The  winds  were  not  strong  enough  to  advect  the  clouds  over  the  land.  The 
cloud  deck  prevented  the  collection  of  satellite  sea  surface  temperature  data  for  the 
region  of  interest  on  the  experiment  day,  September  8,  1989. 

E.  INSTRUMENTATION  CHARACTERISTICS 
AND  PERFORMANCES 

The  relevant  performance  of  assets  which  aided  the  collection  of  data  in  the 
simultaneous  experiments  are  described  in  the  first  part  of  this  section.  The  proce¬ 
dures  of  data  acquisition  and  initial  sorting  are  described  in  the  second  part  of  this 
chapter. 
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The  sensors  which  were  used  in  the  Marina  experiment  on  September  8, 1989 
include  remote  and  in  situ  assets  measuring  synoptic  spatial  variations  and  Eulerian 
and  Lagrangian  temporal  variations  of  the  ocean  surface,  respectively. 

1.  Platforms 

•  NASA/JPL  DC-8  airborne  INSAR  and  SAR  imaging  radar  remotely  sensing 
ocean  surface  velocity  fields  and  surface  intensity  reflectivity. 

•  R/V  Ricketts  for  deploying  and  tracking  Lagrangian  floating  drifters  following 
the  very  near  surface  currents.  High  accuracy  Mini-Ranger  navigation  system 
(~  5m  RMS  error)  for  positioning  the  vessel,  and  hence  drifter  positions, 
was  installed  on  board  of  R/V  Ricketts  in  addition  to  the  inherent  Loran-C 
navigation  system  with  lower  accuracy  (~  65  m  RMS  error). 

•  Light  aircraft  Piper  Cherokee  Arrow  (Figure  6.15)  for  detecting  and  observ¬ 
ing  the  wooden  sheet  floating  drifters,  which  were  undetectable  initially  by 
conventional  means  (radar,  visible)  from  the  surface  platform  (R/V  Ricketts). 

2.  In  situ  Assets 

•  Ocean  wave  directional  spectra  and  surface  wind  data  acquired  from: 

-  Shallow  water  pressure  gauge  array  located  at  a  depth  of  16  m  off  Marina 
Beach  at: 

36°42'  North  Latitude 
121°48.9'  West  Longitude 

-  Anemometer  on  Marina  Beach  located  at  36°42'  N.  Lat,  121°48.9A  W. 
Long  from  which  directivity  of  the  Bragg  waves  may  be  inferred. 

•  Lagrangian  drifters  following  the  very  near  surface  currents  tracked  by  the 
Motorola  Mini-Ranger  navigation  system. 
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•  NAVSTAR.  Global  Positioning  System  (GPS)  extremely  accurate  positioning 
system  for  ground  truthing  the  geographic  coordinates  of  SAR  and  INSAR 
images  (see  Chapter  VII.  B3). 


3.  Platform  Characteristics 
a.  NASA/JPL  DC-8 

The  NASA/JPL  DC-8  (Figure  6.13)  has  a  multifrequency,  multipo¬ 
larization  airborne  SAR  system  operating  at  L,  C  and  P-Band  (Held  et  al.,  1988) 
with  a  full  digital  recording  capability  and  an  inertial  navigator  system  for  measur¬ 
ing  three  dimensionally  the  aircraft  attitude.  Detailed  performance  of  the  system 
characteristics  are  given  in  in  Table  6.3. 

At  L  and  C-Band  the  radar  can  operate  in  an  along  track  interfer¬ 
ometer  mode  by  using  the  physically  separated  fore  and  aft  antenr.as  mounted  on 
the  left  looking  side  of  the  aircraft  fuselage.  In  this  experiment,  the  L-Band  inter¬ 
ferometer  was  used  in  vertical  polarization  since  the  ocean  reflectivity  at  moderate 
incidence  angles  20° -70°  (Bragg  resonant  backscatter)  is  more  intense  as  compared 
with  horizontal  polarization. 

The  receiver  has  no  Sensitivity  Time  Control  (STC)  or  Automatic 
Gain  Control  (AGC)  capabilities.  Instead  a  wide  dynamic  gain  range  (30-56  dB) 
is  implemented  to  compensate  for  signal  saturation.  This  is  achieved  using  8  bit 
Analog  to  Digital  Converters  (ADC)  which  are  clocked  together  to  ensure  registered 
samples.  Due  to  the  high  digitizing  rate  of  the  ADC’s  (45  MHz),  real  data  samples 
(the  imaginary  part  is  not  recorded)  are  produced  at  a  total  data  rate  of  20  to  60 
MBytes/second.  The  data  from  each  receiving  channel  is  multiplexed  together  to 
form  an  extended  “range  line”  for  recording  onto  any  of  the  three  available  High 
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Density  Digital  Recorders  (HDDR)  at  a  rate  of  80  MBits/s  (see  Fig.  6.15).  No  high 
resolution  processing  capabilities  exist  onboard  the  aircraft  (Held  et  al.,  1988). 

b.  R/V  Ricketts 

R/V  Ricketts  (see  Figure  6.16)  is  a  35  foot  light  craft  suited  for  day 
cruises  on  Monterey  Bay.  The  vessel  is  owned  and  operated  by  Moss  Landing  Marine 
Laboratories.  For  her  size  the  R/V  Ricketts  is  reasonably  equipped  with  naviga¬ 
tion  and  deck  handling  equipment  to  function  as  a  field  station  for  oceanographic 
observations. 

4.  In  Situ  Assets  Characteristics 

a.  Shallow  Water  Wave  Array 

The  Marina  shallow  water  wave  array  consists  of  four  pressure  gauges 
configured  in  a  square,  six  meters  on  a  side,  aligned  6°  off  true  North  (see  Figure 
6.16).  The  pressure  sensors  are  located  0.5-1  m  above  a  mean  bottom  depth  of 
16.2  m.  The  northern  sensor  is  a  parascientific  model  2100AS,  digiquartz  pres¬ 
sure  sensor  which  produces  output  voltage  signal  proportional  to  the  absolute  water 
surface  elevation.  The  other  three  pressure  sensors  are  based  on  a  piezo-resistive 
gauge  transducer  manufactured  by  Kulite  Semiconductor.  The  array  was  developed 
to  provide  wave  directivity  in  shallow  water.  The  azimuthal  resolution  of  the  ar¬ 
ray  depends  on  processing  methods,  and  was  4-6  degrees  for  the  range  of  resolved 
wave  periods  between  5-20  sec.  Data  from  the  shore  station  is  acquired  and  ini¬ 
tially  screened  by  a  computer  based  control  station  located  at  Scripps  Institution  of 
Oceanography,  La  Jolla,  California. 

b.  Wind  Speed  and  Direction  from  Skyvane  Wind  Sensor 

A  2102  Skyvane  four  bladed,  low  threshold  propeller  wind  sensor  is 
mounted  on  the  beach  at  Marina  at  an  altitude  of  27  m  from  MSL  (see  Figure  6.17). 
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TABLE  6.3:  NASA/JPL  DC-8  Airborne  Imaging  Radar  Parameters 
Employed  During  Marina  Experiment,  September  8,  1989 


Data 

L-Band 

Frequency,  MHz 

1237-1260 

Center  frequency,  MHz 

1248.75 

Wavelength,  cm 

24 

Chirp  length,  fxs 

11.25 

PRF/V,  Hz  /  m/s 

1.32 

Band  width,  MHz 

20 

Peak  power,  kW 

6 

Azimuthal  beamwidth,  deg 

8 

Vertical  beamwidth,  deg 

44 

Antenna  beamcenter  gain,  dB 

18.3’ 

Azimuthal  pixel  spacing,  m 

3.03  x  4 

Slant  range  spacing,  m 

6.6 

Polarization 

VV 

Receiver  dynamic  range  (dB) 

30-56 

Number  of  bits  per  sample 

8 

Number  of  looks 

1 

Available  number  of  pixels 
(range  x  aximuth) 

750  x  1024 

FWD  Art 


RADAR  CONTROL 
OAIN,  CALIBRATION 
MOOR,  FLIGHT  PLAN 


Figure  6.14:  NASA/JPL  airborne  radar  and  recording  block  diagram. 
The  L-band  interferometric  mode  was  used  for  the  acquisition  of  data  in 
this  study  (from.  Held  et  al.,  1988). 


Figure  6.15:  Light  aircraft  Piper  Cherokee  Arrow  participated  in  the 
Marina  experiment,  September  8,  1989. 


103 


Figure  6.16:  R/V  Ricketts  participated  on  the  final  experiment  on 
September  8,  1989. 

The  accuracy  of  wind  speed  measurement  is  ±  0.5  m/sec  for  wind  speeds  less  than 
15  m/sec.  The  accuracy  for  wind  direction  is  ±  1%  of  the  full  scale  (±180°). 
c.  Lagrangian  Drifters 

Two  types  of  drifters  were  tested  in  our  experiments  to  measure  near 
surface  currents.  In  the  preliminary  experiment  on  July  24,  1989,  a  passive  radar 
free  floating  drogue  was  used.  It  was  tracked  by  a  surface  navigation  radar  installed 
on  R/V  Point  Sur  and  positioned  accurately  by  the  Mini-Ranger  navigation  system. 
A  prototype  test  of  various  drogue  types  were  tested  on  June  1,  1989.  For  further 
details  see  Appendix  C. 

Based  on  the  results  of  the  preliminary  experiments,  a  Lagrangian 
drifter,  made  of  a  2  cm  thick  plywood  sheet  and  0.4  x  1.2  m  in  dimension  was  used 
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for  ground  truthing  in  the  analysis  of  this  study.  The  dimensions  of  the  drifters 
were  selected  after  considering  the  measurements  by  Stewart  and  Joy  (1974)  and 
Goldstein  et  al.  (1989),  in  which  radar  measurements  of  near  surface  currents  were 
compared  with  simultaneous  measurement  by  drifters.  At  L-band  the  corresponding 
ocean  resonant  Bragg  waves  (21-25  cm  at  30°  viewing  angle)  mostly  feel  currents 
in  the  upper  2  cm  of  the  ocean  surface.  Therefore,  the  drifter’s  thickness  insured 
measurement  of  the  near  surface  flows  experienced  by  the  Bragg  backscattering 
waves  while  the  horizontal  dimension  reduced  the  sensitivity  to  vertical  motion 
of  short  waves  (Goldstein  et  al.,  1989).  The  drifters  were  painted  red  with  an 
identification  number  on  both  sides  in  case  of  rolling  over  in  the  surf  zone  (see 
Figure  6.19). 

d.  Mini-Ranger  Navigation  System 

Vessel  positions,  and  hence  drifter  positions,  were  determined  using 
the  Motorola  Mini-Ranger  Falcon  484  system.  The  system  is  widely  used  for  precise 
positioning  in  a  variety  of  nearshore  marine  applications.  The  system  consists  of  a 
transceiving  master  station  borne  on  a  vessel,  and  up  to  four  stationary  reference 
stations  (see  Table  6.4)  along  the  shore.  The  range  to  each  reference  station  is  com¬ 
puted  at  the  master  station  based  on  the  measured  round  trip  travel  time  of  the 
microwave  transmitted  pulse.  The  system  operates  on  the  principle  of  a  pulse  coded 
transponder.  The  master  and  reference  stations  transmit  on  different  frequencies  in 
the  microwave  band,  5510  and  5480  MHz  to  avoid  interference.  The  standard  system 
has  a  maximum  range  of  37  km  (20  nm)  with  reference  station  gain  antenna  of  13 
dB,  and  a  master  station  omni-directional  antenna  gain  of  6  dB.  Ranges  from  up  to 
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Figure  6.17:  Marina  shallow  water  wave  array.  Topi  hardware  configu¬ 
ration;  Bottom :  geometric  orientation. 


Figure  6.19  Plywood  drifter  used  to  Lagrangianly  follow  surface  cur¬ 
rents  during  the  Marina  experiment,  September  8,  1989.  Top :  perspec¬ 
tive  view;  Bottom  apparent  drifter. 


TABLE  6.4:  Station  and  Equipment  Data 


Horizontal  Datum: 


North  American  Datum  of  1983  (NAD  83) 


Computational  Grid:  Universal  Transverse  Mercator  (UTM) 

Zone  10,  central  meridian  123°  W 

Grid  Scale  Factor 

Used  in  Computations:  0.999743 


Source  of  Station 
Positions:  • 


Surveys  by  NPS  Hydrographic  Sciences 
1988-1989 


Reference  Station  Positions 
(see  Figure  6.6): 


Nick- 

Name 

N  Latitude 

W  Longitude 
(deg  min  sec) 

UTM  Coordinates 
Easting  (m) 
Northing  (m) 

Elevation 
MSL  (m) 

WATS 

36  54  26.S95 

602947.7 

23.1 

121  50  39.819 

4085231.4 

HAYS 

36  38  33.826 

607621.3 

137.2 

121  47  45.895 

4055915.3 

DOP 

36  36  5.601 

600434.0 

10.0 

RM1 

121  52  37.479 

4051260.0 

MBARI 

36  37  6.868 

597861.5 

30.7 

RM 

121  54  20.156 

4053118.5 
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four  reference  stations,  with  an  rrns  triangulation  position  of  ~  5  m,  were  recorded 
onboard  (see  next  section). 

F.  DATA  ACQUISITION  AND  INITIAL  DIGITAL 
PROCESSING 

In  this  section,  the  initial  processing  methods  of  the  data  acquired  in  the 
Marina  experiment  on  September  8,  1989  are  discussed.  The  processing  of  the 
airborne  radar  data  at  JPL,  the  shallow  water  wave  array  data  initially  processed  at 
Scripps  Institute  of  Oceanography,  and  the  processing  of  the  drifter  and  GPS  data 
acquisition  are  described. 

1.  Airborne  Radar  Data 

The  onboard  high  density  digitally  recorded  radar  data  is  processed  into 
a  two  dimensional  imaging  map  at  JPL.  The  Airborne  Imaging  Radar  (AIR)  ofF-line 
processor  is  run  on  a  VAX  11/730  minicomputer  assisted  by  an  array  processor  FPS- 
5210.  The  AIR  processor  is  able  to  provide  two  main  products:  a  survey  product 
and  a  high  resolution  product.  The  survey  processor  is  a  quick  look  processor, 
scanning  the  raw  data  tape  to  determine  the  quality  of  the  data  and  whether  the 
data  is  suitable  for  futher  processing.  The  high  resolution  processor  provides  a  high 
quality  complex  image.  The  functions  which  were  used  to  process  the  data  in  the 
present  study  are  as  follows  (Nguyen,  1987,  Held  et  ah,  1988). 
a.  Range  Compression 

The  range  complex  matched  filter  to  the  transmitted  signal  is  gen¬ 
erated  from  a  real  range  offset  signal  whose  spectrum  is  the  complex  conjugate  of 
the  transmitted  signal.  In  order  to  suppress  sidelobes,  the  range  reference  function 
is  weighted  by  a  cosine  function.  The  range  compression  is  done  by  inverse  Fourier 
transforming  the  product  of  two  Fourier  transform  sequences  of  the  match  filter  and 
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the  data.  In  order  to  reduce  spectrum  overlap,  the  data  were  subjected  to  a  four-fold 
decimation  by  applying  a  bandpass  prefilter  in  the  azimuth  dimension  centered  at 
the  Doppler  centroid  frequency.  This  decimates  the  azimuth  resolution  by  a  factor 
of  4  relative  to  the  nominal  SAR  resolution. 

b.  Azimuth  Transform 

The  temporal  range  history,  that  is,  the  phase  history  of  a  station¬ 
ary  synthesized  target,  is  a  function  of  the  instantaneous  range  difference  between 
the  airborne  radar  and  the  target.  Therefore,  the  range  of  this  target  can  be  approx¬ 
imated  as  a  quadratic  function  of  time  over  the  synthetic  radar  aperture.  This  effect 
is  known  as  range  curvature,  or  range  migration,  and  has  to  be  corrected  so  that  the 
energy  of  a  point  target  will  be  contained  only  in  range  bias  parallel  to  the  azimuth 
dimension.  This  correction  is  done  in  the  frequency  domain.  The  range  compressed 
data  are  first  read  into  azimuth  direction  format  and  then  Fourier  transformed  using 
a  8192  point  FFT. 

c.  Range  Cell  Migration  (RCM)  Correction 

The  RCM  correction  algorithm  rearranges  the  data  so  that  the 
quadratic  phase  history  effect  in  the  range  dimension  is  mitigated. 

d.  Azimuth  Compression 

The  data  is  compressed  in  the  azimuth  dimension  by  convolving 
the  data  with  an  azimuth  matched  filter.  This  operation  is  done  in  the  frequency 
domain  by  taking  the  FFT  of  the  azimuth  matched  filter,  multiplying  it  with  the 
FFT  of  the  range  compressed  data,  and  taking  the  inverse  FFT  to  produce  the  one 
look  complex  image. 

For  each  of  the  three  legs  processed  for  this  study,  four  separate  complex 
image  files  of  the  same  scene  for  each  antenna  are  received  from  the  high  resolution 
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AIR  processor.  The  size  of  these  files  are  750  range  lines  with  1024  azimuth  pixels 
in  each  line.  The  initially  processed  pixel  spacing  are  approximately  6.7  m  in  the 
viewing  direction  (range)  and  3  m  in  the  azimuth  direction.  In  order  to  reduce 
speckle  noise,  the  corresponding  pixels  of  each  set  of  the  four  files  are  averaged 
incoherently  in  azimuth.  The  resulting  image  from  each  antenna  has  therefore  a 
degraded  azimuth  pixel  size  by  a  factor  of  4,  namely  12  m.  These  two  sets  of 
files  are  further  registered  and  combined  pixel  by  corresponding  pixel  to  obtain 
a  phase  and  intensity  INSAR  image.  The  phase  of  each  pixel  in  the  combined 
interferogram  is  the  phase  difference  between  the  corresponding  pixels  obtained 
by  the  two  antennas,  while  the  magnitude  image  represents  the  product  of  the 
amplitudes  of  the  corresponding  pixels.  A  brief  description  for  obtaining  the  INSAR 
images  is  depicted  in  Appendix  A. 

Further,  the  aircraft  yaw  and  pitch,  which  cause  sideways  motion  errors 
were  subtracted  from  the  interferogram  resulting  in  pixel  phases  which  represent 
spatial  ocean  surface  velocity  fields  only. 

The  image  processed  from  each  flight  leg  (see  Figure  4.5)  was  received 
from  JPL  on  computer  compatible  6250  BPI  magnetic  tapes.  The  analysis  of 
these  images  was  done  at  the  Air-Ocean  Interactive  Digital  Environmental  Analysis 
(IDEA)  Lab,  which  is  a  remote  sensing  image  processing  laboratory  for  geophysical 
data  at  the  Naval  Postgraduate  School. 

2.  Pressure  Array  Data 

Time  series  records  of  2048  samples  of  pressure  data  from  all  four  chan¬ 
nels  (sensors)  were  processed  for  the  required  experimental  period  1000-1500  PDT. 
These  records  were  examined  for  data  gaps  and  signal  interference.  Some  of  the 
relevant  recognized  errors  which  were  taken  care  of  in  the  initial  preprocessing  by 
Scripps  Institute  of  Oceanography  are:  Spikes  —  Any  value  of  the  time  series  that 
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is  more  than  five  times  the  standard  deviation  from  the  mean  is  considered  as  a 
spike.  Spikes  are  replaced  with  the  previous  value  in  the  time  series.  Absence  of 
zero  crossing  —  a  wave  time  series  (after  detranding  and  demeaning)  that  does  not 
cross  the  zero  mean  level  for  a  specified  number  of  points  is  considered  unacceptable. 
Filtering  —  the  tidal  component  is  removed. 

The  initially  processed  data  by  Scripps  Institute  of  Oceanography  was 
recorded  on  a  9  track  magnetic  tape.  Further  data  analysis  was  done  at  the  Naval 
Postgraduate  School  (NPS)  and  described  in  the  next  chapter. 

3.  Lagrangian  Drifter  Data 

Ranges  from  up  to  four  Mini- Ranger  shore  stations  were  logged  onboard 
the  R/V  Ricketts  on  a  mini-computer  at  5-second  intervals  while  the  vessel  ma¬ 
neuvered  slowly  alongside  each  drifter.  The  positioning  records  indicate  that  vessel 
speeds  while  passing  the  drifters  were  typically  less  than  2  m/sec.  Hence,  vessel  po¬ 
sitions  were  measured  at  a  spacing  of  10  meters  or  less  alongtrack  (during  the  Group 
One  deployments,  data  logging  was  at  10-second  intervals  which  implies  positioning 
data  about  every  20  meters  alongtrack). 

Times  when  each  drifter  was  abeam  the  Mini-Ranger  antenna  were  hand- 
logged  for  later  interpolation  of  vessel  position  from  the  Mini-Ranger  data  file.  The 
drifters  always  were  passed  down  the  starboard  side  at  a  distance  of  about  2  to  3 
meters  from  the  vessel  (antenna). 

Drifter  positions  were  assumed  to  be  the  same  as  the  vessel  position, 
linearly  interpolated  between  successive  Mini-Ranger  fixes  to  the  time  when  the 
drifter  was  abeam.  This  procedure  does  introduce  some  uncertainty  in  drifter  po¬ 
sition,  because  (1),  the  drifter  actually  was  2  to  3  meters  off  the  antenna  position, 
and  (2),  the  interpolation  assumes  the  vessel  followed  a  straight  line  during  the  5  to 
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10  seconds  between  recorded  fixes.  However,  these  assumptions  do  not  practically 
degrade  the  surface  velocity  estimates  as  discussed  in  the  next  chapter. 

4.  Global  Positioning  System  (GPS)  Data 

A  set  of  eleven  ground  points,  visible  in  the  SAR  images,  were  selected 
and  positioned  to  determine  their  ground  coordinates  in  latitude  and  longitude. 
Typical  points  were  road  intersections  and  fence  line  corners.  Because  the  pixel 
size  on  the  images  being  used  was  on  the  order  of  201  meters,  it  was  felt  that  the 
positioning  technique  did  not  need  geodetic  survey  quality  precision.  A  single  GPS 
receiver  was  ut>ed  with  the  antenna  set  as  close  as  practical  to  the  imaged  point. 
Once  the  receiver  had  acquired  the  satellites  and  the  readings  stabilized  (about  15 
minutes),  a  single  position  and  PDOP  value  (precision  indicator)  was  recorded  and 
the  equipment  moved  to  the  next  point. 

The  GPS  receiver  was  operated  in  the  point  positioning  mode,  tracking 
at  least  four  satellites.  It  should  be  noted  that  the  observations  were  made  in  early 
December  1989,  using  Block  I  satellites  before  “selective  availability”  was  activated, 
which  otherwise  could  have  degraded  the  accuracy  of  this  point  positioning  technique 
significantly.  The  standard  error  in  ground  coordinates  was  estimated  to  be  30 
meters,  based  on  prior  experience  using  GPS  for  geodetic  surveys  in  the  area. 

The  ground  coordinates  obtained  from  the  GPS  receiver  were  latitudes 
and  longitudes  on  the  World  Geodetic  System  1984  (WGS  84)  datum,  essentially  the 
same  as  the  North  American  Datum  of  1983  (NAD  83).  Knowledge  of  the  datum 
was  important  because  the  bathymetric  data  was  provided  on  an  earlier  datum, 

NAD  27,  and  the  shift  between  the  two  is  significant. 

'The  rectified  image  is  of  size  1024  x  1024  pixel  (see  Chapter  VII.  B)  monitored  cn  a  screen  of 
512  x  512  (lines  x  columns).  Therefore,  the  azimuth  and  range  pixel  size  are  doubled  (e.g.,  12.1 
x  2  m  in  azimuth)  due  to  this  constrain. 
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Elevations  of  the  ground  points  also  were  needed,  to  compute  a  small 
correction  to  the  range  ( y )  pixel  numbers.  The  elevations  were  scaled  from  a  to¬ 
pographic  map  of  the  area,  specifically  the  MARINA  quadrangle  published  by  the 
United  States  Geological  Survey.  Image  coordinates  for  each  ground  control  point 
were  scaled  from  compressed  versions  of  the  three  images.  It  was  necessary  to  use 
compressed  images  in  order  to  fit  the  original  1024  x  1024  pixel  format  on  the 
512  x  512  display  available.  The  standard  error  in  picking  image  coordinates  was 
estimated  to  be  three  pixels  (0.6%)  on  the  512  pixel  format. 

Of  the  eleven  ground  control  points  originally  selected,  five  were  dis¬ 
carded  because  of  unresolved  blunders,  most  probably  mistaken  identification  of  the 
ground  point  in  the  image.  This  left  a  set  of  six  ground  points  to  use  in  the  adjust¬ 
ments.  A  listing  of  these  points  and  their  coordinates  is  given  in  Table  1,  Appendix 
C,  along  with  other  constants  pertinent  to  the  adjustment  process. 

5.  Summary 

This  chapter  has  described  the  methods,  asset  performance  and  data 
collection  of  the  Wave  Field  and  Surface  Currents  Marina  1989  experiment.  The 
digital  remote  sensing  data  sets  of  the  L-band  INSAR  and  SAR  images  from  the 
orthogonal  and/or  opposite  legs  are  further  processed  to  analyze  surface  currents 
and  wavefields.  These  results  will  also  be  compared  with  conventional  data  sets 
processed  from  in  situ  observations  of  drifters  and  wave  array  data  (see  Chapter 
VII). 
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VII.  DATA  ANALYSIS  METHODS  AND 

RESULTS 


A.  INTRODUCTION 

Methods  by  which  the  raw  data  were  processed  to  obtain  radar  remote  sensing 
spatial  images  and  in  situ  temporal  records  were  described  in  the  previous  chapter. 
In  this  chapter,  the  SAR  and  INSAR  digital  images  are  further  processed  to  obtain 
two-dimensional  wave  number  spectra  of  dominant  waves,  nearshore  surface  current 
vectors,  and  estimations  of  the  decorrelation  time  of  the  backscattering  ocean  sur¬ 
face.  INSAR  wavenumber  spectra  (wavelength  and  direction)  are  compared  for  the 
first  time  with  in  situ  wave  array  spectra  and  theoretical  refraction  and  shoaling 
predictions.  The  INSAR  measured  ocean  surface  current  field  is  compared  with 
Lagrangian  drifting  buoys  following  surface  currents.  The  results  are  in  agreement 
and  demonstrate  the  feasibility  of  INSAR  to  measure  ocean  surface  wave  spectra 
and  currents. 

B.  SAR  AND  INSAR  IMAGE  PROCESSING 

The  SAR  initial  processing  at  JPL  basically  generates  the  two-dimensional 
digital  image  with  a  number  of  inherent  geometric  and  radiometric  distortions.  Usu¬ 
ally  a  number  of  corrections  needs  to  be  applied  in  order  to  have  a  correct  represen¬ 
tation  of  the  imaged  scene.  The  major  post  processing  methods  which  were  applied 
to  the  data  presented  in  this  study  are  described  below.  Digital  image  data  sets 
obtained  from  flight  paths  (Figure  6.5)  in  the  Marina  experiment  September  8, 1989 
are  presented. 
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1.  Geometric  Rectification 

Since  the  image  has  equal  pixel  spacing  in  the  slant  range,  the  data 
have  to  be  interpolated  in  order  to  get  equal  pixel  spacing  in  the  ground  range 
format.  The  raw  image  of  750  lines  (in  range  direction)  of  equal  slant  spacing  and 
1024  pixels  in  each  lines  (in  azimuth  direction)  is  first  scanned  to  find  the  number 
of  initial  unprocessed  lines  in  the  image  provided  by  JPL.  Next,  the  meaningful 
nonzero  lines  are  linearly  interpolated  with  equal  slant  range  pixel  spacing  into  1024 
equal  ground  range  spacing. 


A\gj)  =  A(Si)  + 


A(Si+l)  -  A{St) 

Si+i  —  Si¬ 


te  -  St) 


i  =  1,2,.  ..,750 
j  =  1,2,.  ..,1024 


(7.1) 


where  A'(gj)  is  the  new  interpolated  pixel  value,  S',-,  gj  are  the  slant  and  ground 
range  coordinates  of  a  given  pixel  with  a  fixed  azimuth  coordinate  before  and  after 
rectification,  and  A(S,)  is  the  original  pixel  value  centered  at  5,-.  Since  the  images 
cover  relatively  small  regions  (~  12  km  x  7  km)  the  efTect  of  Doppler  shift  due  to 
the  Earth’s  rotation  which  might  skew  the  image  in  the  azimuth  dimension  is  not 
considered. 


2.  Radiometric  Rectification 

As  described  in  the  previous  chapter,  each  image  processed  data  set 
obtained  from  the  three  different  flight  lines  in  the  Marina  experiment  consists  of: 
two  complex  SAR  images  and  two  interferometric  1024  x  750  phase  and  amplitude 
images.  The  amplitude  or  brightness  of  the  INSAR  image  is  similar  in  a  way  to  the 
conventional  SAR  image  except  that  it  was  obtained  by  taking  the  modulus  of  the 
product  of  corresponding  complex  pixel  images  from  the  two  antennas.  The  image 
pixel  brightness  /,•  is  a  function  of  the  reflectivity  or  the  backscatter  cross  section  a 
of  the  corresponding  imaged  surface  area 


h  =  /(«), 


(7.2) 
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where  / (cr,-)  is  a  function  of  the  radar  transmitted  power,  the  antenna  gain  pattern, 
the  effect  of  range  attenuation,  and  the  transfer  function  of  the  radar  receiver  and 
the  post  processor.  For  an  airborne  imaging  radar  (low  altitude,  small  swath)  it  is 
assumed  that  the  main  contribution  to  brightness  attenuation  is  the  range  parame¬ 
ter,  which  attenuates  the  two  way  transmitted  signal  by  the  fourth  power  R~4.  In 
order  to  mitigate  this  effect  in  the  amplitude  (or  intensity)  images,  each  range  line 
is  normalized  in  the  image  pixel  brightness  by 

*"-(£)*  (7'3) 
where  Ri  is  the  slant  range  to  the  first  line  and  /?,-  is  the  range  of  line  i  (t  = 
1,2,. ,  1024).  This  operation  yields  a  more  uniform  brightness  image  in  the  range 
direction.  Figure  7.1  illustrates  two  brightness  images  from  the  northern  flight  path 
360°  before  and  after  radiometric  rectification;  the  difference  between  the  two  images 
is  clearly  observed.  Figure  7.2  shows  a  rectified  brightness  image  from  the  western 
and  eastern  flight  lines,  legs  270°  and  90°. 

3.  Computation  of  Image  Coordinates 

Ground  coordinates  for  each  flight  path  for  SAR,  and  hence  INSAR, 
images  were  found  by  computing  parameters  of  a  best-fit  transformation  between 
image  and  ground  coordinates  using  a  small  set  of  GPS  surveyed  ground  points 
visible  in  each  image.  A  least-squares,  minimum-variance  adjustment  technique 
was  used  to  solve  for  the  five  parameters  describing  a  simple  rescaling  rotation  and 
translation  transformation: 

k'x  =  pixel  size  scale  in  the  azimuth  direction  (along  track); 
k'y  =  pixel  size  scale  in  the  range  direction  (cross  track); 
a  =  geographic  azimuth  of  the  flight  line; 


117 


Eo  =  easting  ground  coordinate  of  image  origin;  and 
No  =  northing  ground  coordinate  of  image  origin. 

Hence,  the  image  coordinates,  expressed  as  (x,y)  pixel  number  pairs,  were  treated 
as  an  orthogonal  coordinate  system  with  constant,  although  different  metric  scales 
in  the  x  and  y  directions,  the  azimuth  and  range  directions,  respectively. 

A  priori  values  for  the  first  three  parameters  k'x)  kty ,  a  were  provided  with 
the  image  data  from  JPL,  or  alternatively  were  computed  from  information  provided. 
These  a  priori  values  were  used  as  additional  “observations”  in  the  adjustment 
process,  and  were  accepted  as  the  best  values  for  computations  using  image  data. 

Consequently,  the  adjustment  process  served  two  purposes.  First,  it 
confirmed  that  the  a  priori  pixel  scale  sizes  and  flight  path  directions  (azimuth 
values)  were  reasonably  correct.  Second,  the  adjustment  estimated  values  for  the 
ground  coordinates  of  the  image  or  origin,  which  was  needed  in  order  to  register  each 
image  against  the  bathymetric  data  set,  surface  current  ground  truth  measurements, 
ray  traces,  and  station  positions  (Figure  7.9). 

4.  Interferometric  Phase  Images 

The  interferometric  phase  images  (interferograms)  are  further  processed 
to  reduce  the  speckle  noise  and  to  obtain  image  pixels  with  equal  size  in  azimuth 
and  range  dimension.  The  phase  image  contains  the  radial  component  of  ocean 
surface  velocity.  This  physical  information  has  to  be  carefully  preserved  throughout 
the  image;  for  instance,  smoothing  methods  of  the  speckle  noise  have  to  consider 
this  requirement.  The  extreme  dark  and  bright  pixels  representing  high  positive 
and  negative  velocities  are  screened  out  carefully  after  determining  the  dynamic 
range  of  image  pixel  values.  The  upper  and  lower  limits  of  this  dynamic  range  is 
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Figure  7.1:  Conventional  synthetic  aperture  radar  image  of  the  Marina 
Beach  region  in  Monterey  Bay;  Top:  before  radiometric  rectification;  Bot¬ 
tom:  after  radiometric  rectification.  The  aircraft  attitude  was  Eastward 
(90°). 
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Figure  7.2:  Conventional  SAR  image  of  similar  area  as  Figure  7.1  with 
aircraft  attitude  in:  Top:  Northward  (360°)  and  Bottom:  Westward  (270°). 
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determined  by  scanning  the  dark  and  bright  pixel  values  over  the  image  relative 
to  the  “normal”  pixel  values.  Any  pixel  image  value  which  is  beyond  or  under  the 
determined  threshold  is  replaced  by  a  mean  pixel  value  of  its  four  or  eight  nearest 
neighbors  in  azimuth  and  range  dimension.  The  interpolation  of  pixels  for  equal 
size  in  both  dimensions  assists  the  extraction  of  wavenumber  in  the  two-dimensional 
wavenumber  spectra  (see  below).  Figures  7.3,  7.4  are  phase  images  corresponding 
to  the  subswaths  illustrated  in  Figure  6.5. 

C.  MEAN  SURFACE  WIND  DIRECTION 
ESTIMATION 

Fortunately,  the  imaged  scene  in  the  Marina  experiment  included  the  Salinas 
River.  The  phase  images  (Figures  7.3,  7.4)  illustrate  the  difference  between  the 
stationary  coast  and  moving  water  surface.  This  contrast  is  significantly  enhanced 
as  compared  with  the  conventional  SAR  image.  The  image  of  the  coast  is  nearly 
uniform  with  the  expected  zero  velocity  except  for  the  moving  surface  of  the  Salinas 
River.  As  can  be  seen  in  Figures  7.3,  7.4,  and  in  more  detail  in  Figure  7.5,  the  river  is 
blocked  off  from  the  ocean  at  this  time  of  the  year  and  constitutes  a  basin  of  standing 
water  in  which  waves  are  only  generated  by  local  wind  with  no  contamination  by 
long  fetch  dominant  waves  or  advecting  currents.  The  apparent  nonzero  surface 
velocity  in  this  basin  was  measured  and  compared  with  Bragg  phase  velocity  and 
wind  direction.  Anemometer  wind  measurements  at  the  time  of  the  radar  airborne 
overflights  taken  close  to  the  shoreline  of  Marina  Beach  at  a  height  of  27  m  show 
clearly  an  eastward  wind  component  of  almost  2  m/sec.  This  indicates  that  for  leg 
360°  the  radar  was  looking  at  upwind  Bragg  waves.  The  measured  mean  surface 
velocity  obtained  from  the  Salinas  River  is  estimated  by  INSAR  phase  difference 
expressed  by 
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Figure  7.3:  Radar  interferograms  off  Marina  Beach  in  Monterey  Bay, 
September  8,  1989.  Brightness  represents  surface  motions.  The  aircraft 
flight  path  is  from  right  to  left  along  the  top  of  the  image  (Northward 
(360°}  attitude).  A  narrow  band  long  crested  swell  from  south-southwest 
togetner  with  a  wider  band  wave  field  from  north-northwest  are  clearly 
observed  as  they  propagate  ashore.  Top :  before  filtering;  Bottom:  after 
Gaussian  filtering. 
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Figure  7.5:  Salinas  River  interferograms:  Top:  for  flight  path  360°;  Bot¬ 
tom for  flight  path  270°. 
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The  radial  velocity  Ur  is  in  this  case  directly  proportional  to  the  net  contribution 
by  the  Bragg  wave  phase  velocity  propagating  towards  and  away  from  the  radar 
in  the  viewing  direction.  The  radial  component  of  the  phase  velocity  contribution 
is  rectified  by  the  local  incidence  angle  to  obtain  surface  components  in  the  range 
direction: 


CP=-r 


Ut 


sm7 


(7.5) 


The  phase  difference  a,  which  can  vary  between  0  and  2n  radians,  is  repre¬ 
sented  by  0-255  image  pixel  gray  shades  with  zero  velocity  at  pixel  value  of  128. 
The  estimated  surface  phase  velocity  Cpsr  in  the  Salinas  River  from  the  image  pixel 
values  is  given  by 

/'PIXVAL  -  128\  AV 

CpSR  =  { - 128 - )  2Bsin7  (7'6) 

where  PIXVAL  is  the  actual  pixel  image  value.  From  linear  theory,  the  phase 
velocity  of  a  surface  gravity  wave  is  taken  to  be: 


(7.7) 


Substituting  the  Bragg  resonant  condition  (Equation  2.6)  into  Equation  (7.6)  results 


in 


Cp  = 


9  X  1 


(7.8) 


2n  2  sin  7 

Using  Equations  7.5,  and  7.7,  the  estimated  apparent  INSAR  Bragg  waves  phase  ve¬ 
locity  is  computed  from  two  observations  and  compared  with  theoretical  prediction. 
The  results  are  shown  in  Table  7.1. 
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TABLE  7.1:  Surface  Bragg  Wave  Phase  Velocity  Estimated  from  INSAR 
Interferogram  and  Linear  Theory  Phase  Velocity 


Flight 

Path 

7  (deg] 

Cpsr  [cm/sec] 

Cpth  [cm/sec] 

360° 

46 

53  ±2 

51 

270° 

33 

— 411  ±  7 

58 

1  The  minus  sign  indicates  Bragg  waves  propagating  away  from  the  aircraft  in  range  direction. 

Based  on  these  results,  the  surface  wind  is  infered  to  have  had  mainly  an 
eastward  component  with  a  weaker  ambiguous  northerly  component.  The  result 
clearly  supports  the  assumption  that  in  flight  path  360°  the  radar  was  looking  at 
an  upwind  surface  roughness  while  in  flight  path  270°  the  radar  looked  at  a  weakly 
crosswind  scene.  This  indicates  that  the  resonant  Bragg  waves  responsible  for  the 
radar  reflectivity  travelled  ambiguously  in  both  directions  (for  flight  path  270°), 
towards  and  away  from  the  radar,  resulting  in  difficulties  in  estimating  their  apparent 
contribution  to  the  observed  radial  velocity  from  the  ocean  surface. 

D.  CORRELATION  OF  OCEAN  SURFACE 

REFLECTIVITY  WITH  WIND  DIRECTION 

The  estimate  of  the  surface  wind  pattern  from  the  Salinas  River  interfero- 
grams  (Table  7.1),  are  further  infered  regarding  the  mean  surface  reflectivity  ob¬ 
tained  from  the  conventional  SAR  amplitude  images  from  the  three  flight  paths. 
The  surface  wind  direction  as  estimated  in  the  previous  section  had  mainly  an  east¬ 
ward  component  with  a  weaker  component  southward.  The  image  intensity,  which 
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is  a  function  of  the  backscatter  cross-section  of  the  corresponding  surface,  should  in 
turn  be  most  intense  when  looking  upwind  and  weaker  when  looking  crosswind. 

Corresponding  mean  image  intensities  from  the  three  flight  paths  were  taken. 
The  results  are  not  surprising  and  are  consistent  with  the  estimates  of  the  Bragg 
phase  velocity  in  Salinas  River. 

2^1  =  4 -4- 7  (7.9) 

1 90° 

=  1 .5  -T-  3  (7.10) 

■*270° 

These  results  clearly  show  that  the  surface  reflectivity  is  strongest  when  looking 
upwind,  weaker  when  the  radar  is  looking  cross-upwind,  and  weakest  when  looking 
cross-downwind.  This  might  be  one  of  the  reasons  for  the  difference  in  the  quality 
of  the  phase  images  of  opposite  flights  (90°  and  270°)  covering  similar  scenes  (see 
Figure  7.3).  The  eastward  90°  interferogram  seems  to  be  less  noisy  because  the 
essential  backscatter  signal  is  sufficient  enough  to  extract  the  phase  information. 
The  radar  signal  obtained  from  the  Bragg  backscattering  waves  is  probably  less  than 
the  minimum  required  in  the  270°  flight  path  to  obtain  a  meaningful  interferogram 
representing  ocean  surface  dominant  waves  and  currents. 

E.  ESTIMATION  OF  SCENE  COHERENCE  TIME 

The  estimation  of  the  radar  scene  coherence  time  of  the  ocean  surface  is 
processed  following  the  method  described  in  Chapter  IV.  D2.  Image  intensities  of 
similar  subscenes  obtained  separately  from  the  two  antennas  and  from  the  product 
of  the  two  were  averaged.  The  scene  coherence  time  was  computed  using  Equation 
(4.61) 

Tc  =  £  (-2  In  AW1'2 
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The  scene  coherence  time  was  computed  for  all  the  three  flight  paths  for  different 
subscenes.  The  corresponding  averaged  subscenes  were  chosen  from  near  to  far 
incidence  angle.  The  results  are  depicted  in  Table  7.2.  For  the  Northern  flight  path 
(360°),  the  viewing  ocean  scene  is  only  between  40°-50°  of  incidence  angle;  the  scene 
coherence  time  first  increases  slightly  (until  46°)  and  then  remains  almost  constant. 
For  further  discussion  of  these  results,  see  Chapter  VIII.  B. 

F.  AMBIENT  SURFACE  CURRENT  ESTIMATES 
FOLLOWING  DRIFTING  BUOYS 

Vessel  positions  were  computed  from  the  Mini-Ranger  ranges  using  a  least- 
squares  technique  called  “variation  of  coordinates”.  Descriptions  of  the  procedure 
can  be  found  in  most  surveying  texts.  The  presentation  in  Cross  (1981)  is  virtually 
identical  to  that  used  in  the  present  analysis. 

Although  the  Mini-Ranger  system  computes  and  logs  a  position  on  line,  all 
range  data  were  reprocessed  in  order  to  eliminate  errors  caused  by  erratic  reception 
from  the  station  at  MBARI  RM  (Table  6.6).  In  fact  all  data  from  that  station  were 
eventually  rejected. 

1.  Positional  Accuracy 

The  two-dimensional  random  errors  in  vessel  position  are  estimated  as 
part  of  the  “variation  of  coordinates”  procedure.  These  errors  are  primarily  func¬ 
tions  of  the  relative  geometry  between  the  ship  and  shore  stations,  and  the  esti¬ 
mated  precision  of  a  single  range  measurement.  This  random  error  is  expressed  as 
a  variance-covariance  matrix  computed  as  follows: 

C  =  o20{AtPA)~'  (7.11) 
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Table  7.2:  Estimates  of  Scene  Coherence  Time  as  a  Function 
of  Incidence  Angle  Obtained  from  Three  Flight  Paths  in  the  Marina 
Experiment,  September  8,  1989 


Mean  Incidence 

Scene  Coherence  Time  r  [msec] 

Angle  7 

leg  90° 

leg  270° 

leg  360° 

25° 

136 

121 

— 

29° 

146 

123 

— 

33° 

139 

134 

— 

36° 

117 

116 

— 

40° 

121 

105 

119 

42° 

114 

97 

134 

45° 

96 

95 

133 

o 

CO 

109 

95 

138 

50° 

112 

94 

140 

Mean 

121 

109 

132.8 

Standard  Deviation,  a 

15.2 

13.2 

7.35 
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where  C  is  the  2x2  symmetric  variance-covariance  matrix  with  the  diagonal  elements 
being  the  variances  of  easting  and  northing  coordinates  and  the  ofT-diagonal  being 
their  covariance.  The  quantity  <Jq  is  the  “a  priori  variance  of  unit  weight”  which 
scales  C.  Practically,  Oq  is  assumed  to  be  unity  in  these  computations  and  checked 
with  an  “a  posteriori''  estimate  derived  from  the  actual  range  measurements  in  each 
fix. 

The  matrix  /l,  is  the  N  x  2  Jacobian  for  the  transformation  from  the 
//-dimensional  range  coordinate  system  to  the  two-dimensional  UTM  coordinate 
system.  In  this  case  N  typically  was  3  for  the  ranges  from  reference  stations  WATS, 
HAYS,  and  DOP  RMl  finally  used  (see  Table  4.6). 

Ay  —  {Xj  -  Sij)/Ri  i  =  1,2,3;  j  =  1,2;  (7.12) 

where  X  is  the  ship  position,  Si  is  the  z-th  station  position,  and  Ri  is  the  horizontal 
distance  on  the  UTM  grid  between  A'  and  Si. 

The  matrix  P  is  the  N  x  N  weight  matrix  which  is  the  inverse  of  the 
variance-covariance  matrix  for  the  N  range  observations.  It  is  assumed  the  range 
measurements  are  all  independent  with  equal  variances  (it  is  assumed  here  tr  =  3 
meter.,  for  any  single  range  measurement).  Hence  P  becomes  a  simple  identity 
matrix  times  the  scalar  quantity  l/o-2,  or  1/9. 

Because  the  drifters  were  offset  from  the  Mini-Ranger  antenna,  an  ad¬ 
ditional  random  error  must  be  included.  The  offset  distance  was  between  2  and  3 
meters,  but  the  direction  was  not  known.  In  theory,  the  course- made- good  derived 
from  successive  fixes  could  be  used  to  estimate  the  ship’s  heading,  and  hence  the 
direction  to  the  drifter.  However,  uncertainty  in  direction  still  would  be  very  large 
because  of  the  relatively  slow  speeds  and  inherent  “crabing”.  It  was  simpler  to  just 
add  a  circular  random  error  of  2  meters  to  the  variance-covariance  matricies  derived 
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above.  It  can  be  shown  this  is  accomplished  by  simply  increasing  the  easting  and 
northing  variances  by  4  m2  each. 

This  additional  circular  random  error  also  serves  to  represent  an  uncer¬ 
tainty  due  to  the  linear  interpolation  assumption.  No  theoretical  approaches  could 
be  found  in  estimating  the  random  error  in  the  interpolation,  however,  based  on 
simulation  runs,  the  maximum  error  due  to  interpolation  between  fixes  of  10  m 
apart  would  be  about  2  m. 

In  Appendix  D,  the  variance-covariance  matrices  computed  for  each  ob¬ 
servation  of  drifters  from  Group  One  and  Two  are  presented.  Only  one  matrix  was 
computed  for  each  group  of  observation  because  the  drifters  were  so  closely  bunched 
that  the  relative  geometry  did  not  change  significantly  and  hence,  the  random  errors 
were  essentially  the  same  for  all. 

In  general,  errors  increased  slightly  closer  towards  the  shore.  Circular 
errors  increased  from  about  11  m  to  13  m  (CE  95%),  going  from  the  most  offshore 
positions  where  Group  One  was  deployed  to  the  nearshore  area  where  Group  Two 
was  recovered.  Also  note  that  the  northing  coordinate  was  determined  more  pre¬ 
cisely  than  the  easting,  by  a  factor  of  about  3.  Both  effects  are  due  to  the  fact  that 
the  reference  shore  stations  lay  essentially  along  a  north-south  line  passing  close  to 
the  experiment  region. 

2.  Accuracy  of  Velocity  Estimates 

Drifter  velocities  were  computed  by  the  difference  of  positions  divided 
by  time  interval.  Standard  error  propagation  techniques  were  used  to  estimate  a 
variance-covariance  matrix  for  the  U  and  V  velocity  components.  The  reader  is 
refered  to  the  MATHCAD  documentation  for  a  description  of  the  math  and  the 
actual  computation  of  the  estimates  given  below. 
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The  error  estimates  of  the  sorted  data  are  small,  of  the  order  of  0.1 
cm/sec  as  shown  in  Table  7.3.  The  errors  depicted  are  small  primarily  because 
of  the  relatively  long  time  intervals  between  drifter  position  observations,  typically 
about  40  minutes.  The  relatively  large  time  period  of  drifter  positioning  is  allowable 
since  the  environmental  conditions  on  September  8, 1989  were  almost  steady  for  this 
time  interval.  Also  note  that  these  error  estimates  assume  that  U  and  V  components 
are  constant  over  the  period  between  drifter  positionings.  In  other  words,  these  are 
estimates  of  the  error  in  the  mean  velocities. 


TABLE  7.3:  Standard  Deviation  of  Velocity  Components  (cm/sec) 


°u 

<JV 

Group  One: 

Obs  1  to  2 

0.13 

0.08 

Obs  2  to  3 

0.19 

0.11 

Obs  3  to  4 

0.22 

0.12 

Group  Two: 

Obs  1  to  2 

0.14 

0.07 

Obs  2  to  3 

0.20 

0.10 

G.  INSAR  MEASUREMENT  AND  GROUND 

TRUTHING  OF  NEARSHORE  OCEAN  SURFACE 
CURRENTS 

INSAR  phase  images  (Figures  7.3,  7.4)  were  further  processed  to  obtain  esti¬ 
mates  of  very  near  ocean  surface  currents.  The  radial  component  of  surface  velocity 
UT  of  each  pixel  is  rectified  from  the  incidence  angle  dependence  resulting  in  the  hor¬ 
izontal  component  of  ocean  surface  velocity  in  the  viewing  (range)  direction.  The 
bias  introduced  by  the  Bragg  wave  phase  velocity  estimated  from  the  Salinas  River 
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(see  Chapter  VII.  C)  is  then  subtracted.  The  resulting  phase  image  is  averaged 
two-dimensionally  to  mitigate  the  quasi-periodic  wave  orbital  velocity.  The  scale  of 
averaging  is  the  orbital  motion  scale  of  the  largest  observed  wave.  All  smaller  wave 
multiples  of  this  scale  are  averaged  out.  The  surface  velocity  image  from  one  flight 
path  represents  only  one  component  of  the  surface  current  vector  in  the  viewing  di¬ 
rection.  Two  orthogonal  flight  paths  observing  the  same  area  are  needed  to  estimate 
the  full  current  vector.  The  phase  image  data  from  flight  paths  360°  and  90°  were 
selected  to  obtain  interferometric  virtual  estimates  of  ocean  surface  current  fields. 
These  estimates  are  compared  with  simultaneous  observations  of  the  Langrangian 
floating  drifters,  and  good  agreement  is  obtained. 

1.  Additional  Considerations  and  Processing 

When  computing  current  estimates  from  INSAR  phase  images,  the  fol¬ 
lowing  additional  considerations  were  taken: 

•  When  using  the  interferometric  SAR  technique  to  measure  surface  velocities, 
one  has  to  consider  sideways  motion  of  the  antennas,  introducing  a  phase  shift 
making  the  imaged  scene  appear  to  move  when  it  is  stationary.  The  prominent 
errors  in  the  image  phase  data  were  corrected  at  JPL  (see  Chapter  VI.  F)  using 
the  Inertial  Navigation  System  (INS)  data  and  land  imagery  as  a  stationary 
reference  (Goldstein  et  al.,  1989). 

•  Mini-Ranger  positions  of  the  floating  drifters  were  transformed  to  image  pixel 
coordinates  using  the  GPS  adjustment  method  (see  section  B3).  Image  sub- 
areas  from  flight  paths  90°  and  360°  containing  the  observed  path  of  the  La- 
grangian  drifters  were  selected  to  match  the  current  estimates  as  spatially 
similar  in  nature  as  possible  to  the  in  situ  observation. 
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•  The  image  constructed  from  the  360°  flight  path  looking  upwind  is  less  noisy 
compared  with  images  constructed  from  flight  paths  90°  and  270°  looking 
crosswinds  (see  section  C,  D).  Therefore,  prior  to  the  computation  of  current 
estimates,  the  phase  image  from  the  eastward  flight  path  was  smoothed.  Pixel 
phase  values,  which  deviate  by  more  than  two  standard  deviations  from  the 
mean  phase  value  of  a  selected  subarea,  were  rejected,  resulting  in  a  less  noisy 
phase  image  for  current  computation. 

An  additional  problem  is  the  appearance  of  slicks  (seen  as  dark  regions 
in  the  conventional  SAR  amplitude  image,  Figures  7.1,  7.2)  in  the  selected  area.  In 
these  areas  the  essential  backscattering  amplitude  is  probably  too  weak,  i.e.,  very 
noisy,  to  extract  a  reliable  interferometric  phase  due  to  surface  motion.  Nearby 
regions  with  no  slicks  were  selected  to  compute  current  estimates  for  areas  observed 
by  drifters  3,  4. 

2.  Comparison  of  INSAR  Current  Estimates  with  Lagrangian 

Drifter  Measurements 

The  simultaneous  in  situ  surface  current  estimates  are  shown  in  Figure 
7.6.  Three  groups  of  five  Lagrangian  drifters  each  were  deployed  close  to  the  time  of 
INSAR  aircraft  overflights  (see  Chapter  VI.  C).  Radar  data  was  initially  acquired  30 
minutes  before  the  predicted  time  of  maximum  flood  tide  current  (1320  PDT).  This 
timing  was  predicted  on  observing  the  stronger  flows  to  make  the  current  estimates 
more  reliable  in  an  area  (Marina  Beach)  which  is  known  to  have  unpredictable  weak 
tidal  currents. 

The  Lagrangian  drifters  (Figure  7.6)  experienced  a  dominant  eastward 
component  and  a  much  weaker  northward  component.  Some  estimates  of  surface 
current  derived  from  interferometric  SAR  data  and  Lagrangian  drifters  are  tabulated 
in  Table  7.4.  It  is  shown  that  the  dominant  eastward  components  are  in  better  agree- 
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LONG 

Figure  7.6:  Direct  observation  of  surface  currents  following  Lagrangian 
drifters  near  Marina  Beach  in  Monterey  Bay,  California.  Th>.  arrows 
point  in  the  direction  of  the  current  and  the  length  of  each  arrow  is 
proportional  to  the  observed  t /eed. 

ment  (correlation  0.9)  than  the  weaker  and  hence,  difficult  to  determine,  northward 
component  (correlation  0.5).  The  velocity  magnitudes  are  most  consistent  with  a 
correlation  of  0.94.  The  reasons  for  these  results  are  assumed  to  be  associated  with 
the  physical  difference  between  the  two  types  of  measurements,  Lagrangian  versus 
spatial  synoptic  observation,  and  are  further  discussed  in  Chapter  VIII.  C. 

The  drifter  velocity  estimates  with  almost  nil  standard  deviations  (au  = 
0.18  cm/s,  aT  =0.1  cm/s)  were  discussed  in  the  last  section.  The  INSAR  mea¬ 
surements  had  a  larger  standard  deviation  of  20-30  cm/sec  computed  from  selected 
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TABLE  7.4:  INSAR  Estimates  and  In  Situ  Observation  of  Ocean  Surface 
Currents  (in  cm/sec).  Subscript  d  Denotes  Drifter  and  s  for  INSAR. 


Drifter 

Zonal 

Component 

Meridional 

Component 

Velocity 

Magnitude 

Number 

Vd 

ua  (360°) 

Vd 

(90°) 

ud 

Ua 

1,2 

20 

5 

21 

m 

3,4 

18 

15 

5 

4 

19 

1111 

5 

18 

16 

4 

5 

18 

6,7 

22 

17 

5 

9 

23 

19 

8 

21 

17 

5 

8 

22 

19 

9, 10 

22 

18 

4 

5 

22 

19 

Average 

20 

17 

5 

6 

21 

18 

rms  error  (cm/sec) 

0.5. 

1.7 

0.5 

subareas,  which  were  averaged  to  estimate  INSAR  current  components.  This  large 
velocity  deviation  is  believed  to  be  primarily  associated  with  the  orbital  velocity  of 
the  quasi-periodic  surface  gravity  waves  which  are  averaged  out.  To  get  an  estimate 
of  the  order  of  magnitude  of  this  assumption,  consider  a  a  deep  water  surface  gravity 
wave  propagating  in  the  cross-track  (range)  direction  of  a  moving  radar  platform. 
For  this  geometry,  the  observed  radial  component  is  the  actual  wave  orbital  veloc¬ 
ity.  For  waves  propagating  in  the  along  track  (azimuth)  direction,  only  the  vertical 
particle  velocity  is  observed.  Therefore,  the  largest  observed  orbital  velocity  in  deep 
water  is  given  by 


ii  -,isV 

Vob-j~ 


(7.13) 


where  II  is  the  wave  height.  The  orbital  velocity  of  a  shoaling  wave  increases  with 
decreasing  wavelength;  hence,  the  expected  dominant  contribution  to  the  observed 
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mean  orbital  velocity  is  from  the  shorter  observed  swell  component  (see  Chapter  IV. 
C). 

For  a  wave  with  L  —  97  m,  period  T'  —  9.1  sec  and  a  rms  wave  height 
Ha  =  0.32  m  (the  wave  height  and  period  are  estimated  from  the  in  situ  power 
spectrum,  see  Figure  7.10)  as  observed  in  the  Marina  experiment,  UtoB  =  0.20 
m/sec.  This  result  is  consistent  with  the  computed  standard  deviation  from  the 
estimated  mean  surface  currents.  It  indicates  that  most  of  the  contribution  to  the 
standard  deviation  is  from  the  quasi-periodic  orbital  velocity  of  dominant  waves 
which  are  averaged  out.  Additional  contribution  to  the  standard  deviation  of  the 
mean  surface  current  component  estimates  is  attributed  to  speckle  noise.  Estimates 
of  current  velocity  vectors  for  radar  and  direct  observed  measurements  are  shown 
in  Figure  7.7.  The  agreement  is  felt  to  be  encouraging.  For  further  discussion  on 
these  results,  see  Chapter  VIII.  C. 


50  CM/S 

Figure  7.7:  Detailed  presentation  of  Figure  7.6  together  with  INSAR 
current  estimates  (dashed  arrows). 
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H.  OBSERVATION  OF  SPATIAL  AND  TEMPORAL 
NEARSHORE  OCEAN  WAVE  SPECTRA 

Shoaling  and  refraction  of  waves  can  be  characterized  by  two  dimensional 
spectral  analysis.  The  two-dimensional  wavenumber  (k  =  ( kx ,  ky ))  spectra  obtained 
by  INSAR  imagery  are  tested  against  the  conventional  measurement  of  directional 
spectra  computed  from  time  series  records  obtained  from  the  pressure  wave  array 
in  the  shallow  water  off  Marina  Beach.  In  addition,  the  two  swell  components 
observed  by  the  shallow  water  array  were  used  as  inputs  to  compute  wave  ray  traces 
using  Dobson’s  numerical  model  (see  Chapter  V.  B).  Two-dimensional  INSAR  wave 
number  spectra  were  computed  for  subscenes  along  the  ray  traces.  The  results 
are  consistent  with  the  refraction  model  outputs  and  with  the  shallow  water  array 
directional  spectra. 

1.  INSAR  Two-Dimensional  Wavenumber  Spectra 

Following  the  description  in  Chapter  V.  C,  a  software  package  was  de¬ 
veloped  to  calculate  two-dimensional  wavenumber  spectra  from  INSAR  images. 

The  rectified  phase  images  illustrated  in  Figures  7.3,  7.4  were  further 
processed  by  linear  interpolation  (see  Equation  7.1)  to  equalize  the  azimuthal  pixel 
size  with  the  range  size  resulting  in  square  pixels.  This  operation  was  done  to 
simplify  the  computation  of  the  wave  lengths  and  directions  of  the  peaks  of  the 
wave  number  spectra. 

Due  to  the  complex  shallow  water  region  in  which  the  experiment  was 
conducted,  relatively  small  subscenes  had  to  be  selected  to  fulfill  the  homogeneity 
requirement  for  spectral  analysis.  Trade  off  is  between  satisfying  the  requirement 
of  homogeneity  and  obtaining  confidence  in  the  spectral  estimates.  To  obtain  con¬ 
fidence  in  the  spectral  estimates,  a  minimum  sea  surface  area  equivalent  in  size  to 
approximately  6-10  wavelengths  on  a  side  must  be  used  to  calculate  INSAR  image 
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spectra  to  provide  reasonable  estimates  of  wavelength  and  direction.  Therefore,  a 
flexible  2D  FFT  code  was  developed  which  allows  a  selection  of  any  number  of  pix¬ 
els  in  either  dimension  resulting  in  rectangle  sub-areas  for  spectral  analysis  of  any 
size.  If  the  number  of  pixels  in  the  sub-area  is  not  a  power  of  2,  the  sub-area  is 
zero  paded  to  the  nearest  power  of  2  number  of  pixels.  This  increases  artificially 
the  effective  resolution  for  wavelength  and  angle  computations.  The  edges  of  the 
subscenes  were  cosine  tapered  to  reduce  spectral  leakage  and  demeaned.  Figure  7.8 
shows  some  of  the  sub-areas  that  were  subjected  to  spectral  analysis.  The  INSAR 
image  spectra  of  the  corresponding  sub-areas  are  illustrated  in  Figures  7.12  to  7.15. 
They  are  smoothed  with  a  Gaussian  smoothing  kernel  providing  spectral  estimates 
with  25-30  degrees  of  freedom. 

2.  Shallow  Water  Array  Directional  Spectra 

Computation  of  shallow  water  array  directional  spectra  was  done  us¬ 
ing  the  Grauzinis  (1989)  method  (see  Chapter  V.  C).  The  available  Fourier  coeffi¬ 
cients  of  the  wave  field  angular  distribution  function  were  computed  using  Equations 
(5.31 )— (5.35).  The  model  distribution  function  that  was  used  to  match  the  Fourier 
coefficients  is  the  wrapped  Cauchy  distribution 

D{0)  =  —  ( - ^ - - 

2n  (  1  —  2m\  cos  0  +  mj 

This  distribution  gives  suitable  narrow  estimates  and  allows  for  bimodal  directional 
spectra. 

The  spectral  analysis  was  performed  for  a  one  hour  time  series  1200- 
1308  PDT  on  September  8,  1989  with  subsequent  FFT’s  computed  over  128  second 
intervals  with  a  sampling  rate  of  0.5  Hz.  The  averaged  spectral  estimates  have  ap¬ 
proximately  64  degrees  of  freedom.  As  mentioned  in  Chapter  VI.  F2,  the  data  were 
detrended  and  demeaned  before  being  Fourier  transformed  to  eliminate  contamina- 
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Figure  7.8:  Radar  interferogram  of  Marina  Beach  area,  September  8, 
1989.  The  indicated  sub-areas  were  subjected  to  spectral  analysis  with 
results  appearing  in  Figures  7.12-7.15. 

tion  of  the  wave  frequency  bands  due  to  leakage  of  low  frequency  signals.  The  wave 
pressure  power  spectrum  is  calculated  and  transformed  into  surface  elevation  and 
surface  slopes  spectra  by  applying  linear  wave  theory  transfer  functions  as  described 
in  Chapter  V.  C,  and  subsequently  converted  from  a  frequency  angle  presentation  to 
a  directional  wavenumber  presentation  via  intermediate  water  linear  wave  equations. 

3.  Comparison  of  INSAR  2D  Wavenumber  Spectra  with  Temporal 
Wave  Array  Spectra 

Wave  refraction  rays  for  the  two  dominant  waves  observed  by  the  Ma¬ 
rina  shallow  water  array  were  constructed  using  the  Dobson  model.  The  rays  at 
these  peak  frequencies  (0.063,  0.109  Hz)  were  propagated  offshore.  The  predicted 
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12 1°  4  9.6*  W 


Figure  7.9:  Observed  ray  traces  —  from  INSAR  spectral  analysis  versus 
predicted  wave  orthogonals  for  the  two  swells  measured  by  the  shallow 
water  array  of  Marina  Beach.  The  numbers  are  the  centers  (stations)  of 
the  2D  FFT  sub-areas. 
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Table  7.5:  Actual  INSAR  2D  Wavenumber  Wave  Array  Spectra  Results  and  Predicted  Refrac¬ 
tion  and  Shoaling  Data 
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wavelength  and  direction  along  the  refracted  ray  (about  every  800  m)  were  compared 
with  the  actual  wavelength  and  direction  of  the  two  dimensional  INSAR  image 
wavenumber  spectra  from  the  Northern  (360°)  flight  path  image  (see  Figure  7.8). 
The  analysis  is  summarized  in  Figure  7.9  and  Table  7.5. 

The  effects  of  tides  and  coastal  currents  on  wavelength  and  propagation 
direction  were  ignored  and  considered  negligible  when  constructing  the  ray  traces. 
The  comparison  between  the  predicted  and  INSAR  results  clearly  demonstrate  the 
feasibility  of  INSAR  in  imaging  ocean  wave  spectra  in  a  complex  region  under  the 
given  experimental  environmental  conditions. 

An  inherent  limitation  in  applying  the  refraction  model  to  INSAR’s  data 
is  that  the  bathymetry  used  in  the  model  is  not  as  accurate  as  would  be  desired. 
The  density  of  depth  data  varied  based  on  various  ship  tracks  and  was  available 
only  every  6  seconds. 

The  in  situ  wave  array  spectral  results  are  more  reliable  for  comparison 
with  INSAR’s  measurements.  The  surface  elevation  power  spectrum  obtained  by  the 
wave  sensors  is  illustrated  in  Figure  7.10.  The  two  observed  peaks  in  the  spectrum 
correspond  to  wave  periods  of  Ti  =  15.9  sec  and  Ti  =  9.1  sec.  The  wavelengths 
(Station  1,  Table  7.5)  are  calculated  at  the  appropriate  sensor  depth  (16  m)  using 
the  linear  wave  dispersion  relation  (Equation  5.3).  The  higher  energy  region  around 
the  short  wave  peak  ( T\ )  shows  a  range  of  spectral  average  (0.1  to  0.17  Hz),  while 
the  region  around  the  lower  energy  peak  ( Ti  long  swell)  is  much  narrower;  thus,  the 
spectrum  indicates  the  existence  of  two  swells  with  the  shorter  wave  field  having  a 
broader  spectrum. 

The  directional  wave  spectrum  computed  using  the  exact  Fourier  co¬ 
efficient  representation  method  (Grauzinis  method)  is  presented  in  Figure  7.11. 
The  shorter  wave  field  propagates  from  west-northwest  to  east-southeast,  while 
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the  longer  swell  propagates  from  south-westsouth  to  east-northeast.  The  wave¬ 
length  estimates  from  the  spectral  analysis  of  the  two  data  are  identical  while  the 
directions  are  in  good  agreement  (see  Table  7.5). 


Figure  7.10:  In  situ  surface  elevation  power  spectrum  measured  by  the 
Marina  shallow  water  array  on  September  8,  1989,  1200-1308  PDT,  at  a 
depth  of  16  m. 

In  order  to  reduce  bathymetric  constraints,  two  typical  sizes  of  image 
sub-areas  were  selected  for  computing  the  2D  wavenumber  spectra.  For  the  shorter 
wave  field,  sub-areas  of  ~  1400  x  700  m  were  selected  for  computing  wavelength 
and  direction  depicted  in  Table  7.5,  while  for  the  long  crested  swell  sub-areas  of 
~  1400  x  1400  m  with  overlap  between  the  sub-areas  were  used.  The  variation  of 
bottom  topography  over  the  subscene  changes  the  wavelength  and  direction.  This 
distortion  is  reflected  in  the  2D  wavenumber  spectra.  Table  7.5  illustrates  the  2D 
wavenumber  spectra  of  INSAR  along  the  ray  traces  depicted  in  Figure  7.9. 

The  compu'ed  spectra  (Figures  7.12-7.15)  illustrates  shoaling  and  re¬ 
fraction  of  the  two  dominant  swells  over  a  wide  imaged  region  (~  12  x  6  km).  For 
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region  A,  B  where  the  ocean  is  relatively  deeper,  shoaling  and  refraction  are  less 
pronounced  compared  to  the  southern  regions  C,  D. 

In  spite  of  the  limitations  mentioned  above,  INSAR  image  spectra  veri¬ 
fied  by  wave  array  measurement  and  compared  with  model  predictions  are  encour¬ 
aging,  and  demonstrate  the  feasibility  of  this  unique  technique  to  estimate  ocean 
wave  spectra  even  in  a  complex  occonographic  region.  For  further  discussion,  see 
Chapter  VIII.  D. 
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Dircctionol  Wove  Number  Spectra  Dn*J 
Morino  Sholtow  Woter  Array 
Sept  0,  1909  1200-1308  POT 


Figure  7.11:  In  situ  directional  wave  spectrum  taken  at  Marina  Beach 
for  the  same  time  interval  as  Figure  7.8;  Top :  3D  representation;  Bottom: 
contour  presentation. 
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Figure  7.12:  Top :  Shallow  and  Bottom :  deeper  water  INS AR  image 
wavenumber  spectra  of  sub-area  Al,  A2  indicated  in  Figure  7.8.  The 
circles  represent  wavelengths  as  depicted  in  Figure  7.11. 
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Figure  7.14;  Same  as  Figure  7.13,  but  for  sub-area  Cl,  C2.  The  top 

image  spectra  is  obtained  from  a  sub-area  close  to  the  Marina  shallow 
water  array. 


VIII.  RESULTS  AND  DISCUSSION 


A.  INTRODUCTION 

The  aim  of  the  experimental  part  of  this  study  has  been  to  verify  the  capabil¬ 
ities  of  INSAR  to  image  near  surface  ocean  dynamics.  The  method  was  to  compare 
ocean  wave  spectra  and  surface  current  measurement  obtained  remotely  by  INSAR 
with  conventional  in  situ  observations. 

A  discussion  of  the  processed  results  from  Chapter  VII  is  presented  in  this 
chapter.  It  is  shown  that  the  interferometric  technique  appears  to  offer  considerably 
improved  imagery  of  ocean  surface  waves  as  compared  to  conventional  SAR.  The 
additional  capability  relative  to  SAR,  in  imaging  ocean  surface  currents,  seems  to 
be  reliable  for  oceanographic  studies. 

B.  TEMPORAL  RADAR  COHERENCY  OF  THE 
OCEAN  SURFACE 

The  fluctuations  in  the  short  wave  reflectivity  of  the  ocean  surface  during  a 
given  observation  time  causes  the  temporal  correlation  function  of  its  reflectivity  to 
have  a  width  rc.  The  azimuth  resolution  pax  is  degraded  as  compared  to  that  of  a 
stationary  scene  (see  Equation  3.2)  to 


In  addition  to  the  degradation  in  azimuthal  resolution  due  to  finite  scene 
coherence  time,  the  SAR  (and  hence  INSAR)  azimuthal  resolution  is  smeared  also 
by  dominant  wave  motions.  The  azimuthal  wave  pattern  translation  and  the  range 
component  of  the  wave  acceleration  during  the  integration  time  primarily  degrades 
the  azimuthal  resolution.  For  an  imaged  swell,  the  wave  pattern  translation  smearing 
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can  be  corrected  by  focus  adjustment  of  the  processing  matched  filter  (Raney,  1980) 
and  therefore  will  not  be  discussed  further.  As  was  noted  in  Chapter  IV.  C,  the 
finite  scene  coherence  time  determines  the  effective  integration  time  and  hence  the 
azimuthal  resolution.  The  L-band  radar  scene  coherence  time  estimation  from  image 
intensities  (Chapter  VII.  E,  Table  7.2)  is  of  the  order  of  0(1O2)  msec. 

The  range  component  of  acceleration  of  a  deep  water  gravity  wave  is  given 

by 

a(x,<)  =  2 H  (  exP  [“£"]  sin(k  •  x  —  wt)  cos  7  (8.2) 

where  H  is  the  wave  height,  T'  is  the  wave  period,  and  z  is  the  water  depth.  At  the 
surface  2  =  0.  The  degradation  in  azimuth  resolution  due  to  wave  acceleration  is 
taken  to  be  (Alpers  and  Bruening,  1986) 

Paz  =  -  ^  ar(x,  i)T  =  |  ^  2 H  sin(k  • x-ut)T  cos  7  (8.3) 

.  For  typical  JPL/DC-8  parameters,  wave  height  of  0.7  m,  wavelength  of  100  m, 
effective  integration  time  of  0.1  sec.  and  incidence  angle  of  30°,  paz  =  0.7  m. 
Therefore,  the  degradation  in  azimuth  resolution  due  to  long  wave  acceleration  is 
neglected. 

It  is  important  to  note  once  again  that  the  loss  in  azimuthal  resolution  due 
to  finite  scene  coherence  time  associated  with  the  random  motion  of  scattering  el¬ 
ements  within  the  resolution  cell,  results  in  a  consistent  azimuth  image  smear  and 
not  defocus,  and  is  therefore  not  recoverable  by  focus  adjustment  or  integration 
time  settings  (Lyzenga  and  Schuchman,  1983).  Therefore,  the  estimation  of  scene 
coherence  time  remains  a  key  parameter  in  determining  the  azimuthal  cutoff  wave¬ 
length  of  SAR  and  INSAR  azimuthal  resolution.  If  the  coherence  time  rc  is  much 
smaller  than  the  post  processing  integration  time  T,  the  azimuthal  resolution  can 
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be  approximated  (Lyzenga  and  Shuchman,  1983)  to  be: 


Paz  « 


XR 

2Vrc 


(8.4) 


Based  on  the  results  depicted  in  Table  7.2  and  Equation  (8.4),  the  degradation 
in  azimuth  resolution  due  to  finite  scene  coherence  time  for  the  three  flight  path 
images  in  the  Marina  experiment  is  illustrated  in  Figure  8.1. 


Figure  8.1:  Degradation  in  azimuth  resolution  for  the  NASA/JPL  DC-8 
L-band  SAR  for  the  three  flight  path  images  in  the  Marina  experiment, 
September  8,  1989. 


These  estimates  of  L-band  scene  coherence  times  are  consistent  with  con¬ 
ventional  radar  observations  (100  msec)  during  the  TOWARD  experiment  (Plant 
and  Keller,  1984).  The  theoretical  estimates  of  the  scene  coherence  time  using  the 
Pierson-Moskowitz  spectrum  (Table  4.1)  gives  larger  estimates  (3-4  times)  com¬ 
pared  with  the  more  realistic  direct  observation  estimates.  The  estimates  of  the 
ratio  K/I(o  in  Equation  (4.61)  for  all  three  flight  path  images  were  close  to  0.9  and 
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larger.  This  result  suggests  that  for  the  estimation  of  the  scene  coherence  time  only , 
the  measurement  sensitivity  can  be  increased  by  increasing  the  effective  baseline 
separation  (9.9  m)  between  the  antennas.  Variation  of  the  effective  baseline  sepa¬ 
ration  will  allow  to  investigate  the  distribution  of  the  temporal  correlation  function 
which  was  assumed  to  be  Gaussian  in  the  present  study. 

For  rougher  sea  state,  it  is  expected  that  the  observed  scene  coherence  time 
will  be  shorter.  This  is  because  the  waves  become  steeper  with  higher  orbital  veloc¬ 
ities  resulting  in  a  larger  velocity  spread  of  the  scattering  elements.  When  break¬ 
ing  occurs,  the  backscattering  mechanism  is  specular  in  nature  (not  Bragg).  The 
backscattering  slopes  are  propagating  at  the  dominant  wave  phase  velocity  (Lyzenga 
et  al.,  1983).  The  transformation  of  the  advecting  velocity  of  the  backscattering  sur¬ 
faces,  from  orbital  to  phase  velocity  when  breaking  occurs,  further  reduces  the  scene 
coherence  time.  A  cusped  region  is  the  surf  zone,  the  estimated  scene  coherence  time 
from  this  region  is  30%  less  relative  to  estimates  from  nearby  sea  surface  regions 
with  similar  incidence  angle. 

The  difference  in  azimuth  resolution  between  three  flight  path  images  is  not 
fully  understood  at  the  moment  and  is  probably  determined  by  the  speckle  and 
random  noise  and  look  direction  (upwind  crosswind).  Based  on  the  estimates  of 
azimuth  resolution  (40-70  m)  in  Figure  8.1,  the  largest  unresolvable  azimuth  waves 
are  of  wavelength  between  80-140  m  (increasing  with  incidence  angle).  The  orbital 
velocity  of  these  waves  modulating  the  fundamental  scattering  elements  mainly  de¬ 
termine  the  finite  scene  coherence  time  (see  Chapter  IV.  C). 

The  degradation  in  azimuth  resolution  is  one  of  the  major  reasons  for  the  less 
pronounced  wave  patterns  in  the  INSAR  images  obtained  from  flight  paths  270°  and 
90°.  The  azimuthal  orientation  of  the  dominant  swells  makes  them  less  resolvable 
with  increasing  range  as  the  azimuthal  resolution  decreases  (Figure  8.1  and  Chapter 
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VIII.  D.  2c).  The  wave  patterns  in  the  90°  and  270°  SAR  images  were  more  affected 
by  velocity  bunching  compared  to  the  360°  SAR  image.  This  might  contribu.  j  to 
additional  azimuthal  degradation  of  the  INSAR  image. 

In  shallow  water  shoaling  and  refraction  probably  also  play  a  role  in  de¬ 
termining  the  scene  coherence  time.  Due  to  refraction,  waves  tend  to  propagate 
perpendicular  to  the  bathymetric  contours  (refraction)  which  are  mainly  aligned 
alongshore,  while  shoaling  results  in  increased  wave  amplitude  and  elliptic  particle 
orbits.  The  change  in  the  radar  observed  orbital  velocity  due  to  shoaling  is  believed 
to  be  small.  The  more  pronounced  effect  is  probably  the  refraction  effect  in  which 
wave  crests  tend  to  align  with  the  bathymetric  contours  in  shallow  water  before 
breaking  occurs.  In  deeper  water,  the  wave  length  and  direction  are  less  affected  by 
the  bathymetry  and  the  direction  of  wave  propagation  might  be  affected  by  local 
winds.  In  the  Marina  experiment  the  mild  wind  (<  2  m/sec)  had  no  influence  on 
the  direction  of  propagation  of  the  larger  unresolved  waves  (80-140  m)  which  are 
of  our  interest  since  these  waves  determine  the  azimuth  resolution.  Therefore,  the 
direction  of  the  larger  unresolved  waves  in  deeper  water  was  unknown.  For  flight 
path  360°,  the  refraction  effect  caused  waves  to  be  range  oriented  near  the  shore 
(smaller  incidence  angles)  and  hence  the  full  magnitude  of  the  orbital  velocity  should 
be  observed  by  the  radar,  resulting  in  a  faster  decorrelation  time  of  the  surface.  An 
indication  of  this  effect  can  be  noticed  through  the  results  depicted  in  Table  7.2 
for  the  360°  image.  Additional  evidence  using  the  present  technique  and  in  situ 
observations  is  needed  to  get  a  better  understanding  of  the  determination  of  scene 
coherence  time  in  shallow  water. 

An  important  application  of  the  knowledge  of  the  scene  coherence  time  is  in 
defining  the  optimum  physical  baseline  separation  between  the  two  antennas  which 
determines  the  sensitivity  of  interferometric  measurements  of  ocean  surface  velocity 
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fields.  For  a  mean  estimate  of  0.1  sec.  of  the  scene  coherence  time,  it  is  required 
that  the  platform  travel  time  through  the  effective  baseline  separation  (9.9  m)  will 
be  less  than  50  msec.  Using  the  expression  for  the  time  delay  for  interferometric 
measurement  A t  =  rc/2  and  L-band  NASA/DC-8  parameters  yields  a  maximum 
baseline  separation  of: 

Bmta  =  At  •  2V  =  21.4  m  (8.5) 

This  result  shows  that  the  actual  baseline  separation  (19.8  m)  used  in  the  Marina 
experiment  was  almost  optimal  for  the  given  environmental  conditions.  For  rougher 
sea  state  and/or  higher  radar  frequencies,  the  scene  coherence  time  will  be  shorter 
and  hence  the  baseline  separation  between  the  antennas  should  be  reduced.  The 
requirement  for  At  =  rc/2  probably  could  be  more  flexible  allowing  At  <  rc  and 
hence  B^  could  be  increased  resulting  in  an  improved  sensitivity  in  the  INSAR 
phase  measurement  (given  rc  >  0.1  sec.).  The  effects  of  increasing  of  the  antenna 
baseline  separation  should  be  inferred  in  future  experiments. 

C.  INTERFEROMETRIC  SAR  REMOTE  SENSING  OF 
OCEAN  SURFACE  CURRENTS 

A  primary  objective  of  the  present  study  is  the  verification  of  the  feasibility 
of  INSAR  to  estimate  reliably  ocean  surface  currents  with  high  resolution  over  large 
regions.  The  Marina  experiment  demonstrates  this  feasibility  by  comparing  INSAR 
estimates  with  almost  simultaneous  in  situ  Lagrangian  drifter  estimates  (see  Chapter 
VII.  G).  A  scatter  diagram  of  current  magnitude  estimates  derived  from  the  two 
sensors  is  illustrated  in  Figure  8.2.  Simple  statistical  analysis  shows  a  significant 
degree  of  linear  correlation  (R  =  0.94),  with  a  small  consistent  bias  of  2.5  cm/sec. 
However,  due  to  the  relatively  small  number  of  compared  estimates  the  statistical 
significance  is  rather  low. 


156 


0  12  15  10  21  24  27  30 


DRIFTER  CURRENT  ESTIMATE  (CM/S) 

Figure  8.2:  Lagrangian  drifter  versus  INSAR  surface  currents  estimates. 

1.  Limitations  of  Measurement  Methods 

In  the  eastward  (90°)  flight  path  image  (from  which  the  meridional  cur¬ 
rent  component  was  estimated)  the  radar  was  looking  crosswinds.  The  estimation  of 
the  ambiguous  phase  velocity  associated  with  the  Bragg  waves  was  derived  from  the 
western  (270°)  flight  path  image  (see  Chapter  VII.  C,  G),  which  included  the  Sali¬ 
nas  River.  Unfortunately,  the  less  noisy  90°  phase  image  did  not  include  a  nearby 
calm  region  like  the  Salinas  River  to  estimate  the  phase  velocity  bias.  The  270° 
image  was  not  selected  for  computation  of  meridional  current  estimates  because  it 
was  the  most  noisy  image  obtained  in  the  Marina  experiment  (see  Chapter  VII.  D). 
Therefore,  the  substruction  of  the  velocity  bias  introduced  by  the  Bragg  wave  phase 
velocity  for  the  meridional  component  was  not  straightforward,  resulting  in  a  less 
accurate  estimation  of  this  component. 
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Another  difficulty  that  was  addressed  in  estimating  current  components 
from  INSAR  phase  images  was  the  presence  of  slicks,  which  are  clearly  seen  in  the 
amplitude  SAR  images  (see  Figures  7.1, 7.2).  These  slicks  are  believed  to  be  observed 
due  to  local  dampening  (straining)  of  the  wind  generated  surface  ripples  which  are 
responsible  for  the  essential  Bragg  radar  backscatter.  Ocean  surface  straining  can 
be  introduced  by  various  oceanographic  phenomenas,  one  probable  driving  force  to 
create  slicks  is  that  the  surface  signature  is  related  to  the  convergence  of  surface 
particle  velocities  associated  with  internal  waves.  The  internal  waves  might  be 
generated  by  tidal  flow  variabilities  over  topographic  features  causing  fluctuations 
in  the  thermocline.  The  local  convergence  zones  introduce  concentration  of  surface 
organics  which  dampen  the  surface  ripples  due  to  increased  surface  tension.  Another 
mechanism  creating  regions  of  surface  convergence  are  Langmuir  cells  in  which  the 
slicks  are  oriented  in  the  direction  of  the  wind.  Because  some  of  the  slicks  are  aligned 
along  the  bathymetric  contours  and  the  wind  was  weak,  it  is  believed  that  they  are 
related  to  internal  waves  rather  than  Langmuir  cells. 

Mean  temporal  Lagrangian  estimates  are  physically  different  in  character 
from  the  synoptic  instantaneous  spatial  observation  of  the  surface  currents,  which 
could  be  described  as  a  grid  of  Eulerian  snapshots.  The  Lagrangian  drifter  trajectory 
would  not  necessarily  follow  the  spatial  current  pattern  measured  instantaneously 
by  INSAR,  which  is  evident  in  Figure  8.2.  Averaging  of  imaged  subareas  make  the 
two  different  current  estimates  more  comparable.  This  difference  might  be  one  of 
the  reasons  for  the  better  agreement  in  the  computed  current  magnitude  estimates 
rather  than  current  direction  (see  Table  7.4). 

2.  Current  Measurement  Conclusions 

The  results  of  the  experiment  in  estimating  surface  currents  with  INSAR 
shown  in  Table  7.4,  Figures  7.7,  and  8.2,  together  with  the  discussion  above,  suggest 
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a  number  of  conclusions  which  are  basically  consistent  with  those  pointed  out  by 
Goldstein  et  al.,  1989. 

•  The  method  seems  to  be  most  reliable  when  the  radar  observation  is  along  the 
direction  of  the  strongest  current  component  (flight  path  360°). 

•  The  magnitude,  rather  than  the  direction  of  current  estimates,  is  more  con¬ 
sistent  with  conventional  observation,  probably  because  of  the  fundamental 
difference  between  the  two  measurement  techniques. 

•  The  interferometric  technique  seems  to  be  more  accurate  when  the  radar  is 
looking  upwind  (and  probably  also  downwind).  In  these  cases,  the  Bragg 
waves  phase  velocity  responsible  for  the  radar  signal  are  less  ambiguous,  and 
the  returned  amplitude  signal  is  strongest  resulting  in  a  better  signal  over 
noise. 

•  The  estimates  of  the  bias  phase  velocity  (when  looking  upwind)  are  consistent 
with  those  expected  theoretically. 

•  The  interferogram  contains  numerous  additional  features  in  the  ocean  flow 
field  which  were  neglected  in  the  present  study.  Spatial  nearshore  current 
patterns  is  one  phenomena  which  should  be  further  analyzed. 

D.  OCEAN  WAVE  SPECTRA 

The  results  of  INSAR  image  spectral  analysis  tested  against  conventional 
wave  array  spectra  estimates  and  refraction  model  predictions  (Chapter  VII.  H) 
demonstrate  that  interferometric  SAR  under  certain  geometric  orientations  appears 
to  offer  considerably  improved  (relative  to  conventional  SAR)  imaging  of  waves  also 
in  complex  nearshore  regions. 
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1.  Results  Comparison 

The  predicted  (Dobson’s  model)  ray  traces  expressing  the  relationship 
between  the  water  depth  h  and  wavelength  L  and  between  the  depth  and  direction 
of  arrival  a  (Figures  8.3  ,  8.4)  are  such  that  the  refraction  model  becomes  insensitive 
to  changes  in  h  when 

h  >  0.4Io,  (8.6) 

where  Lq  is  the  deep  water  wavelength.  For  the  two  observed  swells  in  the  Marina 
experiment  this  requirement  results  in  a  depth  of  54  m  for  the  shorter  swell  (9.1 
sec)  and  163  m  for  the  longer  swell.  The  western  edges  of  the  Northern  flight  path 
(360°)  INSAR  image  which  was  mainly  subjected  to  spectral  analysis  (Figure  7.8) 
represents  a  mean  depth  of  53  m.  This  indicates  that  both  waves  experience  shoaling 
and  refraction  through  the  entire  imaged  scene. 

The  constructed  wave  refraction  diagram  (Figure  7.9)  is  consistent  with 
the  INSAR  derived  orthogonals.  As  mentioned  previously,  wave  propagation  an¬ 
gles  can  be  sensitive  to  other  oceanographic  phenomena  such  as  tides  and  coastal 
currents.  Wavelengths  are  also  subject  to  these  effects  but  to  a  much  lesser  extent 
(Shuchman  and  Kasischke,  1981,  Hayes,  1980).  These  effects  were  not  considered 
in  our  refraction  model.  The  observed  weak  currents  (Chapter  VII.  G)  and  the 
consistent  INSAR  spectral  results  consistent  with  model  predictions,  indicate  that 
these  effects  were  minor  during  our  experiment. 

Figure  8.5  presents  a  plot  of  wave  directions  c*i,  a2  determined  from 
INSAR  2D  FFT  versus  oi,  a2  calculated  from  Dobson’s  refraction  model.  The 
plot  shows  a  significant  degree  of  linear  correlation  (0.88)  for  the  long  crested  swell 
which  experience  more  significant  refraction  compared  to  the  shorter  observed  swell 
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Figure  8.3:  Relationship  between  water  depth  h  and  predicted  wave¬ 
length  in  Monterey  Bay  based  on  Dobson’s  refraction  model:  Top:  for 
the  shorter  observed  swell  Lq  =  129  m;  Bottom:  for  the  long  crested  swell 
Lq  —  395  m. 
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Figure  8.4:  Relationship  between  water  depth  h  and  predicted  direction 
of  wave  arrival  a  in  Monterey  Bay  based  on  Dobson’s  refraction  model. 
Top:  for  the  shorter  swell  L0  =  129  m$  Bottom :  for  the  longer  swell  L0  = 
395  m. 
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(correlations  of  0.43).  The  reason  for  the  better  agreement  of  the  narrow  band  long 
crested  swell  is  probably  because  this  swell  is  more  monochromatic  in  character 
relative  to  the  shorter  swell  and  hence  more  comparable  to  the  theoretical  model, 
which  is  based  on  the  analysis  of  monochromatic  waves  as  input 

The  results  for  INSAR  ohserved  wavelengths  compared  with  model  pre¬ 
dictions  are  depicted  in  Figure  8.6.  The  linear  correlation  is  better  than  0.9  for  both 
swells.  This  is  mainly  because  shoaling  is  less  sensitive  than  refraction  to  oceano¬ 
graphic  constrains  and  remote  sensing  imaging  geometry.  The  measured  wavelength 
from  INSAR’s  image  spectra  compared  with  the  in  situ  wave  array  wavelength  esti¬ 
mates  and  model  predictions,  ic  in  better  agreement  than  wave  direction  (see  Table 
8.1). 

It  is  felt  that  more  advanced  spectral  analysis  techniques  like  maximum 
entropy  might  improve  the  data  analysis  results  in  shallow  water.  This  and  other 
techniques  offer  the  potential  of  reducing  the  number  of  wavelengths  required  in  a 
sub-area  which  is  subjected  to  spectral  analysis  (Schuchman  and  Kasischke,  1981). 
Together  with  more  accurate  bathymetric  information  the  potential  of  improving 
the  analysis  of  refraction  and  shoaling  in  shallow  water  can  be  obtained. 

TABLE  8.1:  Deviation  of  INSAR  Image  Spectra  Wavelength  and  Direc¬ 
tion  from  Wave  Array  and  Model  Prediction  Output 


Wave  Direction 

Wavelength 

T  =  9.1  sec. 

T  =  15.9  sec. 

T  =  9.1  sec. 

T  =  15.9  sec. 

Correlation 

0.43 

0.88 

0.96 

0.99 

rms  error 

6.2° 

6.4° 

3.9  m 

4.7  m 
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Figure  8.5:  Predicted  refraction  versus  INSAR  2D  wavenumber  spectra 
angles  for  the  two  observed  swells  (9.1  sec.,  15.9  sec.). 
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Figure  8.6:  Predicted  shoaling  versus  INSAR  2D  wavenumber  spectra 
angles  for  the  wave  array  in  situ  observed  swells  (9.1  sec.,  15.9  sec.). 
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2.  Effects  of  Airborne  SAR  System  on  the  Study  of  Ocean  Surface 
Waves 

Since  the  INSAR  technique  is  based  on  SAR  imagery,  some  key  fun¬ 
damental  shortcomings  of  SAR  in  imaging  ocean  surface  waves  were  described  in 
Chapter  III.  D.  Through  examples  from  the  acquired  remote  sensing  data  the  effect 
of  these  shortcomings  on  INSAR  wave  imagery  is  discussed  in  this  section.  The 
known  SAR  coherent  motion  effects  include  velocity  bunching,  acceleration  defo¬ 
cus  and  scene  coherence  time  limitation.  The  noncoherent  scene  motion  effects, 
which  are  only  important  for  relatively  slow  scanning  platforms  like  aircraft,  include 
scanning  distortion  and  smearing  (look  misregistration  for  multi-look  processing), 
a.  Scanning  Distortion 

Scanning  distortion  arises  due  to  the  relative  motion  between  the 
projected  translation  of  the  mean  pattern  of  image  reflectivity  by  the  wave  phase 
velocity  C,  and  the  imaging  radar  platform  velocity  V.  Therefore,  this  relative 
velocity  effects  during  the  observation  time  is  in  general  proportional  to  the  ratio 
C/V. 

The  image  patterns  should  be  corrected  if  the  projected  azimuthal 
component  of  the  phase  velocity  is  significant  compared  to  the  platform  velocity.  In 
this  case  the  image  contrast  can  be  improved  by  adjusting  the  processor  focus  by 
C'/2.  The  spectral  energies  have  to  be  corrected  in  the  azimuthal  direction  by  kC/V. 
Using  this  method  suggested  by  Vachon  et  al.  (1988),  the  inherent  180°  ambiguity 
(for  all  imaging  radars)  in  the  wave  peak  directions  can  be  resolved.  Artifact  spectral 
contours  obtained  from  opposing  flight  directions  will  diverge  after  this  correction, 
while  the  actual  contours  will  converge.  Using  Equation  (3.16)  and  the  computed 
INSAR  wavelength  and  direction  from  spectral  analysis  (see  Table  7.5),  the  scanning 
distortion  effect  is  estimated.  Longer  waves  with  higher  phase  velocity  are  relatively 
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more  distorted.  Therefore,  the  longest  swell  (310  m)  computed  from  station  4  (360° 
image,  Table  7.5)  was  selected  to  demonstrate  the  scanning  distortion  effect.  The 
calculated  wavelength  distortion  is  of  the  order  of  2%  and  5°  in  direction.  For 
shorter  waves  observed  in  360°  image  where  the  waves  were  mainly  range  oriented, 
the  scanning  distortion  for  both, swells  is  practically  negligible.  Figure  8.7  illustrates 
INSAR  image  spectra  computed  from  sub-areas  in  the  wave  array  and  in  deeper 
water  centered  2000  m  west  from  the  wave  array  zone  for  the  eastern  90°  flight 
path.  Shoaling  and  refraction  effects  are  observed  for  the  shorter  swell  only  (peak 
PI)  since  the  long  crested  swell  is  less,  much  less  energetic  in  deeper  water  and  was 
not  observed.  In  shallow  water  the  long  crested  swell  becomes  more  energetic  due 
to  shoaling  and  a  relatively  weak  peak  with  a  wavelength  of  212  m  and  direction  of 
74°,  14°  norther  from  the  computed  ^direction  260°)  wave  array  directional  spectra, 
was  measured.  Due  to  the  azimuthal  orientation  of  the  swells  in  the  90°  flight  path 
the  wavelengths  are  distorted  by  almost  6-7%,  which  is  relatively  a  large  effect 
compared  with  the  corresponding  360°  flight  path  image  in  which  these  swells  had 
a  range  orientation.  Correction  for  scanning  distortion  will  make  these  swells  more 
comparable  to  the  wave  array  wavelength  estimated  (192  m,  Table  7.5). 

In  deeper  water  a  relatively  large  rotation  of  the  shorter  swell  to¬ 
wards  range  direction  is  observed.  This  additional  rotation  (also  of  the  long  swell 
in  shallow  water)  is  not  well  understood  at  the  moment.  The  noncoherent  scene 
motion  effects  are  negligible  for  spacecraft  since  the  ratio  C/V  is  small, 
b.  Velocity  Bunching 

For  azimuthaly  oriented  waves,  the  velocity  bunching  imaging  mech¬ 
anism  is  significant.  Though  the  INSAR  image  formation  is  not  based  on  this 
mechanism  it  is  formed  from  complex  SAR  images,  which  might  be  distorted  due  to 
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Figure  8.7:  INSAR  Eastern  90°  flight  path  image  spectra  computed  from: 
Top :  Shallow  water  sub-area  close  to  the  Marina  wave  array;  Bottom: 
Deeper  sub-area  centered  2000  m  west  from  the  shallow  sub-area. 
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bunching  effects.  For  INSAR  imagery  it  is  felt  that  velocity  bunching  might  con¬ 
taminate  the  image  formation.  This  imaging  mechanism  becomes  less  significant  as 
the  wave  direction  become  more  range  oriented,  the  wavelength  increases  and/or 
the  wave  height  decreases.  Based  on  the  theory  described  in  Chapter  III.  D  and 
INSAR  orbital  velocity  estimates,  the  velocity  bunching  effect  is  evaluated. 

The  fractional  amplitude  of  the  azimuthal  image  shift  is  estimated 
using  Equation  3.13.  The  wave  amplitude  is  computed  from  the  INSAR  image 
data,  using  linear  wave  theory  and  simple  geometry  considerations.  The  mean 
estimated  orbital  velocity  obtained  from  360°  INSAR  image  for  the  9.1  sec.  swell 
varies  between  20-23  cm/sec.  This  velocity  is  consistent  with  the  rough  estimate 
of  Uob  using  Equation  (7.13).  The  shorter  swell  was  selected  for  a  test  case  since 
it  has  a  relatively  larger  orbital  velocity  and  hence  experience  a  stronger  bunching 
effect.  The  imaged  area  is  considered  as  intermediate  water  depth,  therefore  the 
orbital  wave  particle  velocity  pattern  is  elliptic  rather  than  circular.  The  equation 
of  this  pattern  is  given  by: 

U 2  W2 

A 72  +  A*tgh*kh  =  1  ^ 


where  A'  =  AgT'/L,  and  U,  and  W  are  the  horizontal  and  vertical  velocity  compo¬ 
nents.  Transforming  Equation  (8.7)  into  polar  coordinates  after  simple  manipulation 
the  wave  amplitude  is  obtained 


A  = 


Lr 

gT' 


(sin2  7 
tgh2kh 


-f  cos2  7 


) 


1/2 


(8.8) 


with 


_ Uob _ 

(sin2  7  sin2  <j>  +  cos2  <f>y /2 


(8.9) 


where  7  is  the  incidence  angle  and  <t>  is  the  angle  between  the  wave  direction  and 
the  radar  platform  attitude.  For  both  swells  estimated  from  INSAR  image  spectra 
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close  to  the  Marina  wave  array,  C^.,  and  C'  are  computed  using  Equations  (8.8), 
(8.9),  (3.13),  and  (3.14).  Table  8.1  shows  that  both  swells  are  far  from  being  crit¬ 
ical.  The  shorter  swell  is  more  influenced  by  velocity  bunching.  The  highest  ratio 
of  C'/C'a.  =  0.3  is  obtained  for  the  eastward  90°  flight  path  due  to  an  increased 
azimuthal  orientation  of  the  imaged  swell.  This  might  be  one  of  the  reasons  for  the 
more  smeared  INSAR  image  obtained  from  the  eastward  flight  path  relative  to  the 
northward  360°  flight  path  image.  The  results  depicted  in  Table  8.1  mainly  show 
that  the  imaged  scattering  elements  redistribution  is  one  to  one  resulting  in  single 
unambiguous  brightness  patterns  in  the  SAR  imagery  representing  the  actual  ocean 
waves. 

TABLE  8.2:  Actual  Velocity  Bunching  Parameters  (0,0^)  Estimates  Ob¬ 
tained  from  INSAR  Imagery  (Northern  Flight  Path  360°) 


Wave 
Period  T' 
(sec) 

Flight 

Path 

7 

<t> 

h 

(m) 

| 

UoB 

(cm/sec) 

A 

(m) 

C'e  r 

C' 

9.1 

360° 

o 

O 

69° 

16 

97 

23 

0.3 

43 

10 

15.9 

360° 

o 

o 

o 

GO 

16 

192 

15 

0.2 

292 

0.4 

9.1 

90° 

45° 

12° 

16 

120 

16 

0.23 

19 

6 

c.  Coherent  Motion  Effect 

The  degradation  in  azimuthal  resolution  due  to  finite  scene  coher¬ 
ence  time  has  a  low  pass  filter  effect  on  the  SAR  image.  It  shifts  the  azimuthal  wave 
number  to  lower  values.  An  example  of  this  effect  is  inferred  by  comparing  the  wave 
array  spectra  (Figure  7.11,  Bottom)  and  the  corresponding  INSAR  image  derived 
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spectra  (Figure  8.7,  Top).  The  in  situ  wave  array  spectrum  clearly  shows  another 
peak  with  a  wavelength  of  almost  70  m  and  direction  of  300°.  For  flight  path  90° 
this  relatively  short  swell  is  basically  azimuthally  oriented  and  is  not  seen  in  the 
INSAR  image  spectra.  It  is  concluded  that  the  wave  patterns  in  the  INSAR  images 
obtained  from  flight  paths  90°  and  270°  are  less  pronounced  mostly  due  to  the  az¬ 
imuthal  smear  caused  by  the  finite  scene  coherence  time  limitation.  An  evidence  of 
this  conclusion  can  also  be  observed  in  Figure  7.4.  At  low  incidence  angles  where 
the  degradation  in  azimuth  resolution  is  relatively  small  waves  are  clearly  observed 
in  the  INSAR  image.  At  higher  incidence  angles  where  the  degradation  is  larger, 
wave  patterns  are  hardly  seen.  This  evidence  illustrates  the  importance  of  the  finite 
scene  coherence  time  in  determining  the  SAR  and  INSAR  capability  to  image  ocean 
surface  waves. 

3.  INSAR  Derived  Spectra  Conclusions 

•  The  experimental  results  show  how  interferometric  SAR  may  be  used  to  image 
ocean  surface  waves.  Given  the  problems  associated  with  imaging  waves  using 
conventional  SAR,  this  unique  technique  represents  a  significant  advance.  The 
agreement  between  INSAR  image  spectra  and  in  situ  derived  spectra  (Figures 
7.11,  7.14)  is  encouraging.  In  spite  of  this  potential,  this  section  illustrated 
some  of  the  problems  that  might  arise  with  an  airborne  INSAR  in  imaging 
ocean  surface  waves. 

•  Noncoherent  scene  motion  effects  are  dependent  upon  CjV .  For  a  relatively 
slow  airborne  SAR  (and  INSAR)  platform  ,  noncoherent  motion  effects  might 
lead  to  wave  imaging  complications  that  include  scanning  smear  and  distor¬ 
tion,  which  might  result  in  ambiguous  image  wave  patterns.  This  effect  is 
neglibile  for  a  fast  moving  platforms  like  spacecraft.  Focusing  adjustment  and 
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scanning  distortion  correction  methods  are  not  reliable  for  multimodal  wave 
fields,  which  is  usually  the  actual  situation  of  the  ocean  surface  (Vachon  et  al., 
1988). 

•  The  finite  scene  coherence  time  degrades  the  azimuth  resolution  which  is  pro¬ 
portional  to  R/V  and  was  shown  to  have  an  important  role  for  ocean  wave 
imaging  for  aircraft  and  spacecraft.  Another  coherent  motion  effect  is  the  ve¬ 
locity  bunching.  The  observed  wave  field  (flight  path  360°)  contained  mainly 
range-oriented  waves.  This  orientation  is  favorable  for  radar  imagery  since  it 
minimizes  image  distortion  due  to  azimuth  image  shift  (velocity  bunching). 
The  imaged  waves  in  the  Marina  experiment  did  not  experience  significant 
coherent  azimuth  image  shift  (Table  8.2). 

•  Several  effects  degrade  the  INSAR  image  wave  patterns  (see  Figure  7.4).  The 
dominant  effects  are  believed  to  be:  the  degradation  in  azimuth  resolution 
due  to  finite  scene  coherence  time,  looking  crosswinds  rather  than  upwind 
and  speckle  and  random  noise.  The  lesser  effects  for  the  Marina  experiment 
conditions  were  scanning  distortion  and  probably  velocity  bunching. 

•  In  shoaling  waters  the  area  of  the  image  being  subjected  to  spectral  analy¬ 
sis  should  be  reduced  to  allow  the  waves  to  be  represented  locally  without 
refraction  and  shoaling  effects. 

E.  SOURCE  OF  SWELLS 

The  Global  Spectral  Ocean  Wave  Model  (GSOWM)  outputs  for  August  31 
and  September  8,  1989  dearly  illustrate  the  possible  sources  and  directions  of  dom¬ 
inant  wave  propagation  (see  Figures  6.8-6.11).  A  GSOWM  2D  spectrum  taken  at 
a  grid  point  about  60  miles  southwest  off  Monterey  Bay  (Figure  8.8)  indicates  a 


171 


9.5  second  northwesterly  wave  system  together  with  a  narrow-band  southerly  swell 
(16  sec.).  The  propagation  of  the  long  crested  southerly  swell  into  Monterey  Bay 
is  consistent  with  the  increasing  importance  of  refractive  features  with  decreasing 
wave  frequencies  in  the  southern  edge  of  Monterey  Bay. 


Figure  8.8:  GSOWM  directional  frequency  spectrum  obtained  from  a  grid 
point  at  35°  N  Lat  and  122.5°  W  Long  on  September  8, 1989.  A  southerly 
16  sec.  swell  and  a  9.5  sec.  northwest  wave  system  are  observed. 


The  characteristics  of  southerly  swells  in  the  northern  hemisphere  are  not 
well  understood.  Munk  et  al.,  in  a  series  of  papers,  established  the  source  as  high 
latitude  storms  in  the  South  Pacific  and  Indian  Oceans.  These  swells  are  usually  less 
energetic  (smaller  amplitude  after  traveling  10,000  km)  and  longer  in  period  than 
swells  originating  in  the  North  Pacific  (Guza,  1986  and  Semdin,  1986).  Based  on  the 
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global  models  of  pressure  fields,  winds  and  spectral  analysis,  it  is  assumed  that  the 
south-southwest  long  crested  swell  observed  in  the  Marina  experiment,  September 
8,  1989,  was  generated  by  a  high  latitude  southern  hemisphere  storm  and  that  the 
west-northwest  shorter  period  swell  was  generated  by  a  northwest  Pacific  storm. 


IX.  CONCLUSIONS 


The  present  study  reveals  the  considerable  ability  of  this  novel  technique 
to  provide  reliable  quantative  spatial  information  about  dynamic  processes  at  the 
ocean  surface.  Realization  of  the  potential  of  Interferometric  SAR  with  spacecraft 
might  make  it  possible  to  monitor  global  ocean  wave  spectra  and  surface  currents. 

The  mean  estimate  of  the  L-band  scene  coherence  time  0  (100  msec)  obtained 
from  INSAR  and  SAR  imagery  during  the  Marina  experiment,  on  September  8, 1989, 
is  consistent  with  the  sparse  known  observational  estimates  obtained  under  similar 
moderate  sea  state  and  ambient  conditions.  This  parameter  is  shown  to  be  the 
dominant  limiting  factor  in  determining  the  SAR  and  INSAR  azimuth  resolution 
and  hence  the  azimuthal  cutoff  wavenumber  for  imaging  ocean  surface  dominant 
waves.  Waves  with  azimuthal  component  shorter  than  0  (100  m)  could  not  be 
observed.  The  azimuth  image  shift  (velocity  bunching)  caused  by  the  coherently 
sensed  radial  component  of  the  orbital  velocity  of  the  long  waves  might  disturb  the 
INSAR  phase  measurement.  It  is  shown  that  the  noncoherent  scanning  distortion 
effect  is  important  if  the  projected  phase  velocity  of  the  imaged  wave  is  significant 
relative  to  the  platform  velocity.  This  effect  is  negligible  for  fast  scanning  platforms 
like  spacecraft.  It  is  concluded  that  range  oriented  waves  are  favorable  for  INSAR 
(and  SAR)  where  the  undesired  motion  effects  like  the  temporal  decorrelation  of  the 
ocean  surface  and  velocity  bunching  do  not  distort  the  image.  Also,  range  waves 
project  a  larger  radial  component  of  orbital  velocity  (which  is  sensed  by  INSAR), 
compared  with  azimuthal  waves.  Additional  theoretical  and  experimental  studies 
are  required  to  quantify  and  to  provide  better  understanding  of  these  issues  in  terms 
of  INSAR  performance  to  image  ocean  surface  wave  fields. 
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The  results  of  this  experiment  in  estimating  ocean  surface  currents,  together 
with  earlier  work  by  Goldstein  et  al.  (1989),  and  for  the  first  time  ocean  wavenumber 
spectra  suggest  some  additional  conclusions: 

•  The  method  is  shown  to  be  most  reliable  when  the  radar  observation  is  along 
the  direction  of  the  strongest  current  component. 

•  The  interferometric  technique  appears  to  be  more  accurate  when  the  radar  is 
looking  upwind  (and  probably  downwind).  In  these  cases,  the  phase  velocity 
bias  from  the  Bragg  waves  responsible  for  the  radar  signal  is  less  ambiguous, 
and  the  returned  amplitude  signal  is  strongest. 

•  Shoaling  and  refraction  effects  are  clearly  observed  in  the  INSAR  phase  images. 
In  shoaling  waters,  the  area  of  the  image  being  interpreted  via  spectral  analysis 
should  be  reduced  to  insure  homogeneous  waves. 

•  The  interferograms  contain  numerous  additional  ocean  surface  features  which 
were  not  pursued  in  the  present  study,  such  as  convergence,  divergence,  and 
vorticity  fields  and  internal  and  other  long  waves. 
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APPENDIX  A 

INTERFEROGRAM  PROCESSING 


The  recorded  SAR  data  from  each  antenna  consists  of  a  set  of  four  files.  Each 
file  consists  of  1024  azimuth  by  750  range  pixels.  PmtJ-  and  Qmij  are  the  corre¬ 
sponding  complex  compressed  pixel  values  in  the  imaging  map  as  processed  from 
the  recordred  data  received  from  each  of  the  two  antennas.  The  subscripts  m,i,j 
are  the  file  number,  row,  and  column  of  a  pixel,  respectively.  7tJ-  is  the  complex 
image  value  of  a  pixel  in  the  interferogram.  A,j,  is  the  magnitude  and  phase  of 
Ar 

Processed  Imaging  Files  from  each  Antenna 


|*-I024  -*| 


Iij  is  derived  by  taking  the  complex  conjugate  product  of  the  sum  of  the  four 
corresponding  original  complex  image  pixels  obtained  from  each  antenna. 
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4 


(Al) 

m=l 

Iij  =  Aij  exp[t$y] 

(A.2) 

The  interferogram  file  as  received  from  JPL  is  processed  to  a  2048  azimuth 
by  750  range  pixels.  The  first  1024  azimuth  pixels  are  magnitudes  of  the  complex 
reflectivity  (i4y)  and  the  other  1024  pixels  are  phases  which  consists  of  information 
of  surface  velocity. 

The  pixel  averaging  in  Equation  A.  1  is  in  expanse  of  azimuth  resolution.  The 
pixel  size  is  the  resulting  interferogram  is  12.1  m  is  azimuth  and  6.6  m  in  the  viewing 
direction. 
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APPENDIX  B 

ESTIMATION  OF  RANDOM  ERROR  IN 
L AGRAN GIAN  DRIFTER  VELOCITY 


A.  THEORETICAL  BACKGROUND 

Drifter  velocity  is  a  simple  function  of  distance  traveled  divided  by  time, 
where  distance  is  found  by  differencing  positions.  Thus,  random  errors  in  the  posi¬ 
tions  do  propogate  into  the  velocity. 

The  velocity,  V,  is  a  function  of  the  starting  and  ending  position,  Xi  and  X2, 
and  the  time  difference  A T  =  Tj  -  7\ 


V  =  F(X1,X2,AT)  (5.1) 

where  V,  Xj  and  X2  are  all  vector  quantities  with  east  and  north  components  given 
by: 

V  =  (u,v)  Xi  =  (*i,yi)  X2  =  (x2,y2)  (5.2) 


Specifically, 

u  =  (x2-  x^/AT  v  =  (s/2  -yi)/AT  (5.3) 

The  error  propagation  equation  can  be  set  up  with  E(V)  being  the  2x2 
variance-covariance  matrix  of  the  velocity  components 


S(V)  =  GS(X1,X2)  gt  = 


(5.4) 


where  G  is  the  2x4  Jacobian  of  partial  derivatives  and  is  defined  as 


G, 


dFi{u,v ) 


i  =  l,2 


dSjiX^X-)  j  =  1,2, 3, 4 


or  simply  in  our  case 


G  = 


10-1  0 
0  1  0-1 


1 

AT 


(5.5) 


(5.6) 
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and  Ex  is  the  4x4  combination  of  variance- covariance  matrix  for  the  variable 
positions  and  is  taken  to  be 


E(XlfX») 


< 

0 

0  ' 

**2 

< 

0 

0 

0 

0 

O’siyi 

0 

0 

crx,v, 

<  . 

(B.  7) 


B.  OBSERVATION  ERROR  ESTIMATES 

Using  Equations  (B.l)  through  (B.7),  the  velocity  error  estimates  are  com¬ 
puted: 

1.  Group  I:  Between  observation  2  and  1 


T\  =  11,49,30  (hours,  minutes,  seconds)  PDT 
T2  =  13,15,44  PDT 


AT 

S(Xa) 

£(Xa) 

S(V) 


86.2  min 

20.07  2.69 
2.69  8.24 

22.76  2.69 
2.69  8.03 


m 


m 


1.6 -lO"2  2  •  10-3 
2 • 10-3  6 • 10-3 


2  -2 
cm  sec 


<ju  —  0.13  cm/sec 


<t„  =  8  ■  10-2  cm/sec 


2.  Group  I:  Between  observations  3  and  2 

T3  =  14, 16, 12  PDT 
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T2  =  13,15,44  PDT 
AT 
S(Xi) 


60.5  min 
25.59  2.7 
2.7  7.88 


E(Xa) 

S(V) 


22.76  2.69 
2.69  8.03 


m* 


m 


3.7  •  10"3  4  •  10~3 

4-10"3  1.2 -IQ"3 


cr<i  — 


0.19  cm/sec 
av  =  0.11  cm/sec 

3.  Group  II:  Between  observations  2  and  1 

Ti  =  12,05,28  PDT 
T2  =  13,36,08  PDT 


AT  =  90.7  min 

S(Xx)  = 

'  27.2  2.13  ' 

2.13  7.67 

S(X2)  = 

'  32.26  1.49  ' 

1.49  7.44 

E(V)  = 

‘  2  •  10“2  10 
10"3  5  •  1 

m 

„2 

m 

-3 

n-3 


2  —2 
cm  sec 


<ru  =  0.14  cm/sec 
(7V  =  7.10-2  cm/sec 

4.  Group  II:  Between  observations  3  and  2 

T3  =  14,42,39  PDT 


2  -2 


cm  sec 


180 


T2  =  13,36,08  PDT 


=  6.65  min 


S(X0  = 


S(X2)  = 


£(V)  = 


33.78  0.24 
0.24  7.33 

32.26  1.49 
1.49  7.44 


4.1  •  10- 


9  •  10- 


=  0.2  cm/sec 


2  -2 
cm  sec 


crv  =  0.1  cm/sec 
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APPENDIX  C 


IMAGE  COORDINATES  GROUND 
TRUTHING  USING  GPS  DATA 

A.  IMAGE  ADJUSTMENT  PREPARATION 

A  modified  Transverse  Mercator  projection  was  used  to  convert  the  latitude 
and  longitude  ground  coordinates  to  a  plane  coordinate  grid  of  eastings  and  nor¬ 
things,  (E,N)  coordinate  system  (Floyd,  1985).  The  projection  was  centered  in 
the  image  area  in  order  to  minimize  grid  distortions.  GRS  80  ellipsoid  values  were 
used  which  are  essentially  the  same  as  for  the  WGS  84  geometric  ellipsoid.  The 
projection  parameters  are: 


Ellipsoid  Semi-Major  Axis, 
an  Inverse  Flattening,  1/f 

6378137.0  meters 
298.257222101 

Central  meridian  (longitude) 
False  latitude 

False  easting 

Scale  factor,  central  meridian 

W  121°  47.0’ 

N  36°  40.0’ 

3000  meters 
1.00000 

The  range  pixel  number  ( y )  of  the  image  control  points  needed  a  small  cor¬ 
rection  for  elevation  distortion  (Ulaby,  1982).  The  corrected  range  pixels,  yc,  were 
computed  as  follows: 

Vc  =  UK  (> Jky(y  “  !)  +  ro)2  +  2#z  -  z1  -  +  1  (C.l) 

where 
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TABLE  C.l:  Ground  and  Image  Coordinate  Data 


Type 

Latitude  (d  m  s)1 
Longitude  (d  m  s) 
Elevation  (m) 

■ 

Plane 

Coordinates 

(m) 

Easting 

Northing 

Image 

Coordinates  (pixels) 

X 

y 

A90  A270  A360 

Bridge 

N  36  43  58.34 

3204.3 

313 

197 

abuttment 

W  121  46  51.77 
6.1 

7347.1 

139 

69 

Road 

N  36  43  15.46 

2630.4 

146 

295 

245 

intersection 

W  121  47  14.89 
15.2 

6025.0 

484 

222 

104 

Fence 

N  36  42  47.87 

4169.5 

82 

360 

corner 

W  121  46  12.88 
30.5 

5174.9 

434 

280 

Fence 

N  36  42  37.74 

3645.5 

104 

339 

296 

corner 

Fence 

corner 

W  121  46  33.99 
33.5 

N  36  42  12.89 
W  121  46  16.03 
42.7 

4862.4 

4091.3 

4096.4 

415 

86 

364 

305 

24 

Road 

N  36  41  15.95 

1568.8 

190 

250 

397 

intersection 

W  121  47  57.65 
3.0 

2341.3 

261 

493 
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IMAGE 

PARAMETERS 

Azimuth  Scale,  kx  | 

[meters/pixel) 

24.20 

24.20 

24.20 

(a  priori  values) 

Range  Scale,  ky  (meters/pixel) 

16.14 

14.54 

12.34 

Flight  azimuth,  i 

a  (degrees) 

90° 

270° 

000° 

IMAGE 

CONSTANTS 

Flying  height,  H  (meters) 

8285 

8273 

8178 

Ground  range  to  1st  pixel,  r0 
from  nadir  (meters) 

2601 

2918 

4015 
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y  =  range  pixel  number  (uncorrected) 
ky  —  a  priori  range  scale  (meters/pixel) 
r0  ~  ground  range  from  nadir  to  y  =  1  pixel  (meters) 
h  «  flight  altitude  (meters) 
z  =  elevation  of  ground  point  (meters) 

The  magnitude  of  these  corrections  was  small,  generally  less  than  2  pixels.  However, 
for  a  few  points  at  elevations  of  30  meters  or  more  near  the  flight  line,  the  corrections 
ranged  from  3  to  5  pixels.  The  correction  always  increased  the  pixel  number,  moving 
the  point  further  from  the  flight  line.  Only  the  control  points  needed  this  correction, 
the  imaged  water  surface  assumed  to  .be  at  zero  elevation  were  not  subjected  to 
elevation  distortion. 

B.  IMAGE  ADJUSTMENT  COMPUTATION 

The  adjustment  can  be  viewed  as  solving  for  the  best  mathematical  mapping 
between  the  image  and  ground  coordinate  systems,  between  (x,y)  and  (E,  N),  sub¬ 
ject  to  assumptions  about  the  underlying  form  of  the  mathematical  transformation. 

The  adjustment  used  a  least-squares  technique  in  which  the  observed  data  are 
weighted  inversely  proportional  to  their  variances,  as  described  in  most  surveying 
or  photogrammetry  texts.  This  weighting  strategy  results  in  a  minimum  variance 
solution,  which  is  generally  considered  the  best  solution  when  using  data  of  unequal 
precision.  The  specific  form  applied  here  is  called  the  “combined  observation  and 
condition  equation”  method  (Uotila,  1986). 

The  assumed  transformation  was  a  simple  scaling,  rotation,  and  translation 
between  the  (x,y)  and  (E,N)  coordinates  systems,  with  the  addition  of  separate 
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scales  in  the  a:  and  y  directions  (i.e.,  the  pixels  were  assumed  to  be  rectangular,  not 
square).  The  assumed  form  of  the  transformation  equations  was: 


Ei  —  E0  +  X{kx  sin  a  +  y{kv  cos  a  —  0  (C.2) 

Nt  -  N0  +  X{kx  cos  a  -  sin  a  =  0  (C.3) 


where 


Ei,  Ni  =  ground  coordinates  of  ith  control  point, 

Xi,  yi  =  image  coordinates  of  tth  control  point, 
kx ,  ky  =  scales  (meters/pixel)  in  x  and  y  directions, 
a  =  azimuth  of  the  flight  line  and  rotation  between 
coordinate  systems(a  =  0  points  to  north), 

E0,  N0  =  ground  coordinates  of  image  origin,  (x,y)  =  (0,0). 


Note  that  the  ar-axis  points  are  opposite  to  the  flight  line  direction  because  of  the 
way  the  image  processing  hardware  loaded  the  image  data.  This  accounted  for  the 
unusual  appearance  of  the  signs  in  these  equations  (C.2,  C.3).  Also  note  that  a  is 
measured  clockwise  from  north,  which  causes  the  role  of  the  sine  and  cosine  factors 
to  be  reversed  in  comparison  with  some  other  presentations. 

In  the  adjustment  terminology,  the  quantities  kx,ky,a,E0,  and  N0  are  pa¬ 
rameters  to  be  solved  for.  Hence  there  were  5  parameters  in  each  image.  For  each 
control  point,  one  pair  of  condition  equations  (C.2  and  C.3)  were  formed.  As  long 
as  the  number  of  condition  equations  exceeded  the  number  of  parameters,  the  least 
squares  solution  was  possible. 

Three  addition  condition  equations  were  used,  reflecting  the  a  priori  estimates 
of  the  parameters  kx,  ky,  and  a: 


kx  kxo  —  0 


(C.4) 
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(C.5) 

(C.6) 


ky  ky0  -  0 

a  -  a0  =  0 

where  kxo,  J fcyo,  and  a0  were  the  a  priori  values  of  the  respective  parameters. 

Five  control  points  were  used  on  images  A90  and  A270,  and  four  control 
points  on  image  A360.  Hence,  13  condition  equations  were  used  for  A090  and  R270, 
and  11  for  A360;  more  than  enough  for  the  least  squares  adjustment.  Weighting 
required  estimation  of  the  standard  errors  (precision)  on  each  observation.  The 
estimated  values  were  as  follows: 


(te  =  =  30  meters 

ox  —  cv  =  3  pixels 
cra0  =  1  degree 
<r*JO  =  0.1  meters/pixel 
°kyo  —  0.2  meters/pixel 

The  adjust  proceeded  as  follows.  Define  the  set  of  condition  equations  as  F(La,  Xa)  = 
0,  where  the  vector  La  represents  the  adjusted  values  of  the  observables  Ni,  a 
i/i,  kxo,  kyo,  a0)  and  the  vector  Xa  represents  the  adjusted  values  of  the  parameters 
(kx,  ky,  a,  E0,  N0).  The  set  of  condition  equations  F(La,  Xa)  has  the  form  of 
Equations  (C.2)  through  (C.6). 

With, 

N  =  number  of  control  points  on  the  image 
R  =  number  of  condition  equations 
U  =  number  of  parameters  =  5 


then, 


R  =  2N  +  3 


(C.7) 
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which  is  the  number  of  condition  equations  in  the  set  F(La,Xa).  The  length  of  L 
is  4 N  +  3.  The  degrees  of  freedom  in  the  adjustment  is  R-U, 

dof  =  2jV  +  3-  5  =  2AT-2  =  2(N  -  1)  (C.8) 

Because  the  set  F(L0tX0)  ^  0  when  L0  is  the  raw  observations  and  X0  is  the 
approximate  values  of  the  parameters,  define  a  non-zero  vector  W,  length  R, 

W  =  F{L0iX0). 

The  adjusted  parameters  are  found  by  iterating  the  following  matrix  equation  until 
the  Xa  vector  converges, 

Xa  =  X0  +  dX  (C.9) 

dX  =  -(A'M-'Ay'A'M-'W  (C.10) 

where  matrix  A  =  dF/dX  ( R  rows  by  U  columns),  matrix  M  =  BCB 1  (R  by  i?), 
matrix  B  =  dF/dL  (R  by  X),  and  matrix  C  =  covariance  matrix  of  observables  ( L 
by  L). 

In  this  case,  C  was  a  simple  diagonal  matrix  of  estimated  variances  of  each 
observed  quantity,  and  provides  the  relative  weighting  of  the  observations  in  the 
adjustment.  The  traditional  weight  matrix  is  called  P,  where  P  =  C~l. 

For  analysis  of  the  quality  of  the  adjustment,  residuals  on  the  observed  quan¬ 
tities  were  computed,  as  a  vector 

V  =  La-L0  =  -CBtM-1(AdXa  +  W)  (C.ll) 

and  the  “ variance  of  unit  weight ”  was  computed  as  the  scaler 


”1  =  (V‘PV)/dof. 


(C.12) 


The  adjustment  also  computed  precisions  for  the  derived  values  of  the  parameters, 
expressed  as  a  covariance  matrix  of  the  parameters ,  C*, 

Cx  =  (A*M“lA)~l  (C.13) 

Each  of  these  adjustment  results,  Equations  (C.9)  through  (C.13),  were  com¬ 
puted  for  each  image.  The  results  are  summarized  in  Table  C.2. 

C.  POST-ADJUSTMENT  IMAGE  SHIFT 

Normally,  an  adjustment  process  like  that  described  would  be  sufficient  to 
accurately  align  an  image  with  ground  coordinates.  In  this  case,  however,  the  ad¬ 
justed  results  might  have  been  degraded  by  the  lack  of  available  control  points.  In 
typical  image  adjustment  problems  the  control  points  are  selected  to  span  the  en¬ 
tire  image  area.  Because  most  of  each  image  was  over  water,  the  control  points 
were  concentrated  in  a  relatively  small  land  portion.  Hence,  the  danger  the  results 
might  not  be  reliable  in  the  water  area.  Essentially,  the  danger  was  one  of  using 
extrapolation  rather  than  interpolation. 

Examination  showed  that  the  adjusted  scale  and  azimuth  values  agreed  rea¬ 
sonably  well  with  their  a  priori  values  —  but  not  exactly  (most  within  1  standard 
deviation).  But  because  these  three  parameters  determined  the  overall  size  and  ori¬ 
entation  of  the  image,  it  was  decided  to  accept  the  a  priori  values  rather  than  the 
adjusted  values.  The  a  priori  values  were  derived  from  flight  navigation  data  and 
known  characteristics  of  the  SAR  instrument.  Consequently,  it  was  felt  they  better 
represented  the  true  values  for  the  enthe  image  area. 

This  left  the  adjustment  to  determine  the  location  of  the  images.  The  im¬ 
age  origin  coordinates  ( E0,N0 )  solved  for  in  the  adjustment  reflected  the  adjusted 
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TABLE  C.2:  Image  Adjustment  Results 


kg  ky 

a 

E0 

No 

Image 

Parameter  (m/pix)  (m/pix) 

(deg) 

(m) 

(m) 

A90 

Initial  24.20  16.14 

90.0 

Adjusted  24.20  16.20 

89.9 

6173 

-1864 

a  (0.04)  (0.07) 

(0.3) 

(35) 

(33) 

Post  adjustment  origin  shift 

-8 

5 

Final  origin  coordinates 

6165 

-1859 

Variance  of  unit  weight  0.14 

Latitude  limits  N  36-39-00  to  N  36-43-28 

Longitude  Limits  W  121-44-54  to  W  121-53-12 

A360 

Initial  24.20  12.34 

000.0 

. 

. 

Adjusted  24.23  12.51 

001.3 

4129 

11988 

a  (0.12)  (0.22) 

(0.8) 

(91) 

(60) 

Post  adjustment  origin  shift 

-384 

68 

Final  origin  coordinates 

3745 

12056 

Variance  of  unit  weight  1.76 

Latitude  limits  N  36-39-49  to  N  36-46-30 

Longitude  limits  W  121-46-30  to  W  121-50-44 

TABLE  C.2:  Continued 


Image 

Parameter 

kx 

(m/pix) 

ky 

(m/pix) 

a 

(deg) 

E0 

(m) 

No 

(m) 

A270 

Initial 

24.20 

14.54 

270.0 

_ 

_ 

Adjusted 

24.19 

14.34 

271.0 

-4416 

9437 

a 

(0.09) 

(0.14) 

(0.7) 

(68) 

(105) 

Post  Adjustment  origin  shift  —69  —58 

Final  Origin  coordinates  -4485  9379 

Variance  of  unit  weight  1.39 


Latitude  limits  N  36-41-03  to  N  36-41-03 
Longitude  limits  W  121-43-42  to  W  121-52-01 

rather  than  a  priori  image  size  and  orientation.  The  further  the  origin  point  from  the 
actual  control  points  used,  the  less  precise  its  ground  coordinates  were  determined. 
To  compensate  for  this  effect,  ground  coordinates  of  the  four  image  corners  were 
computed  in  two  ways,  using  (1)  the  adjusted  and  then  (2)  the  a  priori  parameter 
values.  The  origin  coordinates  were  then  shifted  to  align  the  centers  of  the  two 
versions.  In  other  words,  any  systematic  errors  were  forced  to  the  edges  of  the 
image  and  minimized  near  the  center. 

The  amount  of  this  post-adjustment  origin  shift  was  barely  significant  on  im¬ 
ages  A90  and  A270,  in  comparison  with  the  standard  error  of  the  adjusted  coordinate 
values  (see  Table  C.2). 

On  image  A360,  however,  this  shift  was  large,  about  180  meters  westward, 
which  was  twice  the  estimated  standard  error.  No  obvious  explanation  for  this 
discrepancy  was  found,  and  subsequent  comparisons  of  the  predicted  wave  gage  po¬ 
sition  with  the  distance  from  the  shoreline  in  the  image  indicated  that  an  additional 
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200  meter  westward  shift  was  needed;  for  a  total  shift  of  about  380  meters  west  of 
the  adjustment  derived  position.  Whatever  the  error,  it  did  not  appear  to  affect  the 
north-south  direction. 
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